16#ifndef dealii_dof_accessor_h
17#define dealii_dof_accessor_h
30#include <boost/container/small_vector.hpp>
39template <
typename number>
41template <
typename number>
43template <
typename number>
46template <
typename Accessor>
54 namespace DoFCellAccessorImplementation
56 struct Implementation;
59 namespace DoFHandlerImplementation
61 struct Implementation;
64 struct Implementation;
70 namespace DoFHandlerImplementation
72 struct Implementation;
86 namespace DoFAccessorImplementation
104 template <
int structdim,
int dim,
int spacedim>
120 template <
int dim,
int spacedim>
130 struct Implementation;
211template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
232 using BaseClass = typename ::internal::DoFAccessorImplementation::
233 Inheritance<structdim, dimension, space_dimension>::BaseClass;
302 template <
int structdim2,
int dim2,
int spacedim2>
309 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
316 template <
bool level_dof_access2>
353 template <
bool level_dof_access2>
389 typename ::internal::DoFHandlerImplementation::
390 Iterators<dim, spacedim, level_dof_access>::line_iterator
391 line(
const unsigned int i)
const;
398 typename ::internal::DoFHandlerImplementation::
399 Iterators<dim, spacedim, level_dof_access>::quad_iterator
400 quad(
const unsigned int i)
const;
449 const unsigned int fe_index =
460 std::vector<types::global_dof_index> &dof_indices,
461 const unsigned int fe_index =
470 const std::vector<types::global_dof_index> &dof_indices,
496 const unsigned int i,
497 const unsigned int fe_index =
507 const unsigned int vertex,
508 const unsigned int i,
509 const unsigned int fe_index =
541 const unsigned int fe_index =
590 std::set<unsigned int>
611 get_fe(
const unsigned int fe_index)
const;
623 "This accessor object has not been "
624 "associated with any DoFHandler object.");
678 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
686 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
719 const unsigned int fe_index =
724 const unsigned int i,
729 const unsigned int vertex,
730 const unsigned int i,
732 const unsigned int fe_index =
739 template <
int,
int,
int,
bool>
748 friend struct ::internal::DoFHandlerImplementation::Policy::
750 friend struct ::internal::DoFHandlerImplementation::Implementation;
751 friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
752 friend struct ::internal::DoFCellAccessorImplementation::Implementation;
753 friend struct ::internal::DoFAccessorImplementation::Implementation;
765template <
int spacedim,
bool level_dof_access>
824 const unsigned int vertex_index,
850 template <
int structdim2,
int dim2,
int spacedim2>
857 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
909 template <
bool level_dof_access2>
941 typename ::
internal::DoFHandlerImplementation::
942 Iterators<1, spacedim, level_dof_access>::line_iterator
943 line(const
unsigned int i) const;
951 typename ::
internal::DoFHandlerImplementation::
952 Iterators<1, spacedim, level_dof_access>::quad_iterator
953 quad(const
unsigned int i) const;
1000 std::vector<
types::global_dof_index> &dof_indices,
1001 const
unsigned int fe_index =
AccessorData::invalid_fe_index) const;
1012 std::vector<
types::global_dof_index> &dof_indices,
1013 const
unsigned int fe_index =
AccessorData::invalid_fe_index) const;
1033 types::global_dof_index
1035 const
unsigned int vertex,
1036 const
unsigned int i,
1037 const
unsigned int fe_index =
AccessorData::invalid_fe_index) const;
1056 types::global_dof_index
1058 const
unsigned int fe_index =
AccessorData::invalid_fe_index) const;
1120 "This accessor
object has not been "
1164 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
1167 const
DoFAccessor<structdim2, dim2, spacedim2, level_dof_access2> &) const;
1172 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
1175 const
DoFAccessor<structdim2, dim2, spacedim2, level_dof_access2> &) const;
1203 const
unsigned int i,
1204 const
types::global_dof_index index,
1205 const
unsigned int fe_index =
AccessorData::invalid_fe_index) const;
1218 friend struct ::
internal::DoFHandlerImplementation::Policy::
1250template <
int structdim,
int dim,
int spacedim = dim>
1268 const int level = -1,
1269 const int index = -1,
1285 template <
typename OtherAccessor>
1295 const unsigned int fe_index =
1306 const unsigned int fe_index =
1326template <
int dimension_,
int space_dimension_,
bool level_dof_access>
1336 static const unsigned int dim = dimension_;
1341 static const unsigned int spacedim = space_dimension_;
1397 template <
int structdim2,
int dim2,
int spacedim2>
1404 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
1486 periodic_neighbor(
const unsigned int i)
const;
1494 neighbor_or_periodic_neighbor(
const unsigned int i)
const;
1507 boost::container::small_vector<
1536 neighbor_child_on_subface(
const unsigned int face_no,
1537 const unsigned int subface_no)
const;
1546 periodic_neighbor_child_on_subface(
const unsigned int face_no,
1547 const unsigned int subface_no)
const;
1577 template <
class InputVector,
typename number>
1598 template <
class InputVector,
typename ForwardIterator>
1601 ForwardIterator local_values_begin,
1602 ForwardIterator local_values_end)
const;
1624 template <
class InputVector,
typename ForwardIterator>
1628 const InputVector & values,
1629 ForwardIterator local_values_begin,
1630 ForwardIterator local_values_end)
const;
1656 template <
class OutputVector,
typename number>
1659 OutputVector & values)
const;
1692 template <
class InputVector,
typename number>
1694 get_interpolated_dof_values(
1695 const InputVector &values,
1697 const unsigned int fe_index =
1761 template <
class OutputVector,
typename number>
1763 set_dof_values_by_interpolation(
1765 OutputVector & values,
1766 const unsigned int fe_index =
1768 const bool perform_check =
false)
const;
1779 template <
class OutputVector,
typename number>
1781 distribute_local_to_global_by_interpolation(
1783 OutputVector & values,
1784 const unsigned int fe_index =
1801 template <
typename number,
typename OutputVector>
1804 OutputVector & global_destination)
const;
1820 template <
typename ForwardIterator,
typename OutputVector>
1823 ForwardIterator local_source_end,
1824 OutputVector & global_destination)
const;
1839 template <
typename ForwardIterator,
typename OutputVector>
1843 ForwardIterator local_source_begin,
1844 ForwardIterator local_source_end,
1845 OutputVector & global_destination)
const;
1853 template <
typename number,
typename OutputMatrix>
1856 OutputMatrix &global_destination)
const;
1862 template <
typename number,
typename OutputMatrix,
typename OutputVector>
1866 OutputMatrix & global_matrix,
1867 OutputVector & global_vector)
const;
1895 std::vector<types::global_dof_index> &dof_indices)
const;
2030 set_dof_indices(
const std::vector<types::global_dof_index> &dof_indices);
2042 update_cell_dof_indices_cache()
const;
2132 friend struct ::internal::DoFCellAccessorImplementation::Implementation;
2136template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2140 return level_dof_access;
2145template <
int structdim,
int dim,
int spacedim>
2146template <
typename OtherAccessor>
2148 const OtherAccessor &)
2151 ExcMessage(
"You are attempting an illegal conversion between "
2152 "iterator/accessor types. The constructor you call "
2153 "only exists to make certain template constructs "
2154 "easier to write as dimension independent code but "
2155 "the conversion is not valid in the current context."));
2163#include "dof_accessor.templates.h"
DoFAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
DoFAccessor< 0, 1, spacedim, level_dof_access > & operator=(const DoFAccessor< 0, 1, spacedim, level_dof_access > &da)=delete
DoFAccessor< 0, 1, spacedim, level_dof_access > & operator=(DoFAccessor< 0, 1, spacedim, level_dof_access > &&) noexcept=default
DoFAccessor(const Triangulation< 1, spacedim > *, const int=0, const int=0, const DoFHandler< 1, spacedim > *dof_handler=0)
DoFAccessor(DoFAccessor< 0, 1, spacedim, level_dof_access > &&)=default
DoFAccessor(const Triangulation< 1, spacedim > *tria, const typename TriaAccessor< 0, 1, spacedim >::VertexKind vertex_kind, const unsigned int vertex_index, const DoFHandler< 1, spacedim > *dof_handler)
DoFAccessor(const DoFAccessor< 0, 1, spacedim, level_dof_access > &)=default
DoFAccessor(const DoFAccessor< structdim2, dim2, spacedim2, level_dof_access2 > &)
typename::internal::DoFHandlerImplementation::Iterators< dim, spacedim, level_dof_access >::line_iterator line(const unsigned int i) const
friend struct ::internal::DoFHandlerImplementation::Policy::Implementation
static constexpr unsigned int space_dimension
const DoFHandler< dim, spacedim > & get_dof_handler() const
void set_mg_dof_index(const int level, const unsigned int i, const types::global_dof_index index) const
DoFAccessor< structdim, dim, spacedim, level_dof_access > & operator=(const DoFAccessor< structdim, dim, spacedim, level_dof_access > &da)=delete
DoFAccessor(DoFAccessor< structdim, dim, spacedim, level_dof_access > &&)
DoFAccessor(const DoFAccessor< structdim, dim, spacedim, level_dof_access2 > &)
void set_mg_dof_indices(const int level, const std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DoFHandler< dim, spacedim >::invalid_fe_index)
TriaIterator< DoFAccessor< structdim, dim, spacedim, level_dof_access > > child(const unsigned int c) const
const FiniteElement< dim, spacedim > & get_fe(const unsigned int fe_index) const
types::global_dof_index vertex_dof_index(const unsigned int vertex, const unsigned int i, const unsigned int fe_index=DoFHandler< dim, spacedim >::invalid_fe_index) const
static constexpr unsigned int dimension
void set_dof_handler(DoFHandler< dim, spacedim > *dh)
DoFAccessor(const DoFAccessor< structdim2, dim2, spacedim2, level_dof_access2 > &)
types::global_dof_index dof_index(const unsigned int i, const unsigned int fe_index=DoFHandler< dim, spacedim >::invalid_fe_index) const
DoFAccessor< structdim, dim, spacedim, level_dof_access > & operator=(DoFAccessor< structdim, dim, spacedim, level_dof_access > &&)
DoFAccessor(const DoFAccessor< structdim, dim, spacedim, level_dof_access > &)=default
DoFAccessor(const Triangulation< dim, spacedim > *tria, const int level, const int index, const DoFHandler< dim, spacedim > *dof_handler)
bool operator!=(const DoFAccessor< structdim2, dim2, spacedim2, level_dof_access2 > &) const
void get_dof_indices(std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DoFHandler< dim, spacedim >::invalid_fe_index) const
void copy_from(const TriaAccessorBase< structdim, dim, spacedim > &da)
typename::internal::DoFHandlerImplementation::Iterators< dim, spacedim, level_dof_access >::quad_iterator quad(const unsigned int i) const
bool fe_index_is_active(const unsigned int fe_index) const
void copy_from(const DoFAccessor< structdim, dim, spacedim, level_dof_access2 > &a)
void set_dof_index(const unsigned int i, const types::global_dof_index index, const unsigned int fe_index=DoFHandler< dim, spacedim >::invalid_fe_index) const
bool operator==(const DoFAccessor< structdim2, dim2, spacedim2, level_dof_access2 > &) const
types::global_dof_index mg_dof_index(const int level, const unsigned int i) const
types::global_dof_index mg_vertex_dof_index(const int level, const unsigned int vertex, const unsigned int i, const unsigned int fe_index=DoFHandler< dim, spacedim >::invalid_fe_index) const
unsigned int n_active_fe_indices() const
unsigned int nth_active_fe_index(const unsigned int n) const
static bool is_level_cell()
void get_mg_dof_indices(const int level, std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DoFHandler< dim, spacedim >::invalid_fe_index) const
typename ::internal::DoFAccessorImplementation::Inheritance< structdim, dimension, space_dimension >::BaseClass BaseClass
DoFHandler< dim, spacedim > * dof_handler
DoFAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
std::set< unsigned int > get_active_fe_indices() const
void set_mg_vertex_dof_index(const int level, const unsigned int vertex, const unsigned int i, const types::global_dof_index index, const unsigned int fe_index=DoFHandler< dim, spacedim >::invalid_fe_index) const
boost::container::small_vector< TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > >, GeometryInfo< dimension_ >::max_children_per_cell > child_iterators() const
DoFCellAccessor(DoFCellAccessor< dimension_, space_dimension_, level_dof_access > &&)=default
void get_active_or_mg_dof_indices(std::vector< types::global_dof_index > &dof_indices) const
void get_dof_values(const InputVector &values, Vector< number > &local_values) const
face_iterator face(const unsigned int i) const
DoFCellAccessor< dimension_, space_dimension_, level_dof_access > & operator=(const DoFCellAccessor< dimension_, space_dimension_, level_dof_access > &da)=delete
const FiniteElement< dimension_, space_dimension_ > & get_fe() const
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > child(const unsigned int i) const
DoFCellAccessor(const Triangulation< dimension_, space_dimension_ > *tria, const int level, const int index, const AccessorData *local_data)
void set_active_fe_index(const unsigned int i) const
void get_dof_values(const InputVector &values, ForwardIterator local_values_begin, ForwardIterator local_values_end) const
void distribute_local_to_global(const Vector< number > &local_source, OutputVector &global_destination) const
const FiniteElement< dimension_, space_dimension_ > & get_future_fe() const
void get_dof_values(const AffineConstraints< typename InputVector::value_type > &constraints, const InputVector &values, ForwardIterator local_values_begin, ForwardIterator local_values_end) const
unsigned int active_fe_index() const
void distribute_local_to_global(const FullMatrix< number > &local_source, OutputMatrix &global_destination) const
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > neighbor(const unsigned int i) const
DoFCellAccessor(const DoFAccessor< structdim2, dim2, spacedim2, level_dof_access2 > &)
DoFCellAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
void clear_future_fe_index() const
void distribute_local_to_global(const FullMatrix< number > &local_matrix, const Vector< number > &local_vector, OutputMatrix &global_matrix, OutputVector &global_vector) const
void distribute_local_to_global(const AffineConstraints< typename OutputVector::value_type > &constraints, ForwardIterator local_source_begin, ForwardIterator local_source_end, OutputVector &global_destination) const
void set_future_fe_index(const unsigned int i) const
void set_dof_values(const Vector< number > &local_values, OutputVector &values) const
void set_mg_dof_indices(const std::vector< types::global_dof_index > &dof_indices)
boost::container::small_vector< face_iterator, GeometryInfo< dimension_ >::faces_per_cell > face_iterators() const
DoFCellAccessor(const DoFCellAccessor< dimension_, space_dimension_, level_dof_access > &)=default
void get_mg_dof_indices(std::vector< types::global_dof_index > &dof_indices) const
void distribute_local_to_global(ForwardIterator local_source_begin, ForwardIterator local_source_end, OutputVector &global_destination) const
~DoFCellAccessor()=default
unsigned int future_fe_index() const
bool future_fe_index_set() const
DoFCellAccessor< dimension_, space_dimension_, level_dof_access > & operator=(DoFCellAccessor< dimension_, space_dimension_, level_dof_access > &&)=default
void get_dof_indices(std::vector< types::global_dof_index > &dof_indices) const
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > parent() const
DoFInvalidAccessor(const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *local_data=0)
typename InvalidAccessor< structdim, dim, spacedim >::AccessorData AccessorData
typename TriaAccessorBase< structdim, dim, spacedim >::AccessorData AccessorData
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
#define DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
static ::ExceptionBase & ExcVectorNotEmpty()
#define DeclException0(Exception0)
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotActive()
static ::ExceptionBase & ExcMatrixDoesNotMatch()
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcInvalidObject()
static ::ExceptionBase & ExcCantCompareIterators()
static ::ExceptionBase & ExcVectorDoesNotMatch()
static ::ExceptionBase & ExcMessage(std::string arg1)
const ::Triangulation< dim, spacedim > & tria