15#ifndef dealii_dof_accessor_h
16#define dealii_dof_accessor_h
28#include <boost/container/small_vector.hpp>
37template <
typename number>
39template <
typename number>
41template <
typename number>
44template <
typename Accessor>
52 namespace DoFCellAccessorImplementation
54 struct Implementation;
102 template <
int structdim,
int dim,
int spacedim>
118 template <
int dim,
int spacedim>
128 struct Implementation;
209template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
230 using BaseClass = typename ::internal::DoFAccessorImplementation::
231 Inheritance<structdim, dimension, space_dimension>::BaseClass;
300 template <
int structdim2,
int dim2,
int spacedim2>
307 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
314 template <
bool level_dof_access2>
351 template <
bool level_dof_access2>
387 typename ::internal::DoFHandlerImplementation::
388 Iterators<dim, spacedim, level_dof_access>::line_iterator
389 line(
const unsigned int i)
const;
396 typename ::internal::DoFHandlerImplementation::
397 Iterators<dim, spacedim, level_dof_access>::quad_iterator
398 quad(
const unsigned int i)
const;
447 std::vector<types::global_dof_index> &dof_indices,
459 std::vector<types::global_dof_index> &dof_indices,
468 const std::vector<types::global_dof_index> &dof_indices,
494 const unsigned int vertex,
495 const unsigned int i,
506 const unsigned int vertex,
507 const unsigned int i,
587 std::set<types::fe_index>
620 "This accessor object has not been "
621 "associated with any DoFHandler object.");
675 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
683 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
715 const unsigned int i,
721 const unsigned int i,
727 const unsigned int vertex,
728 const unsigned int i,
736 template <
int,
int,
int,
bool>
744 friend struct ::internal::DoFHandlerImplementation::Policy::
746 friend struct ::internal::DoFHandlerImplementation::Implementation;
747 friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
748 friend struct ::internal::DoFCellAccessorImplementation::Implementation;
749 friend struct ::internal::DoFAccessorImplementation::Implementation;
761template <
int spacedim,
bool level_dof_access>
846 template <
int structdim2,
int dim2,
int spacedim2>
853 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
905 template <
bool level_dof_access2>
937 typename ::
internal::DoFHandlerImplementation::
938 Iterators<1, spacedim, level_dof_access>::line_iterator
939 line(const
unsigned int i) const;
947 typename ::
internal::DoFHandlerImplementation::
948 Iterators<1, spacedim, level_dof_access>::quad_iterator
949 quad(const
unsigned int i) const;
996 std::vector<
types::global_dof_index> &dof_indices,
997 const
types::fe_index fe_index =
numbers::invalid_fe_index) const;
1008 std::vector<
types::global_dof_index> &dof_indices,
1009 const
types::fe_index fe_index =
numbers::invalid_fe_index) const;
1029 types::global_dof_index
1031 const
unsigned int vertex,
1032 const
unsigned int i,
1033 const
types::fe_index fe_index =
numbers::invalid_fe_index) const;
1052 types::global_dof_index
1054 const
types::fe_index fe_index =
numbers::invalid_fe_index) const;
1116 "This accessor object has not been "
1117 "associated with any DoFHandler object.");
1160 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
1163 const
DoFAccessor<structdim2, dim2, spacedim2, level_dof_access2> &) const;
1168 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
1171 const
DoFAccessor<structdim2, dim2, spacedim2, level_dof_access2> &) const;
1199 const
unsigned int i,
1201 const
types::fe_index fe_index =
numbers::invalid_fe_index) const;
1213 friend struct ::
internal::DoFHandlerImplementation::Policy::
1245template <
int structdim,
int dim,
int spacedim = dim>
1263 const int level = -1,
1264 const int index = -1,
1280 template <
typename OtherAccessor>
1300 const unsigned int i,
1321template <
int dimension_,
int space_dimension_,
bool level_dof_access>
1331 static const unsigned int dim = dimension_;
1392 template <
int structdim2,
int dim2,
int spacedim2>
1399 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
1502 boost::container::small_vector<
1532 const unsigned int subface_no)
const;
1542 const unsigned int subface_no)
const;
1572 template <
class InputVector,
typename number>
1593 template <
typename Number,
typename ForwardIterator>
1596 ForwardIterator local_values_begin,
1597 ForwardIterator local_values_end)
const;
1619 template <
class InputVector,
typename ForwardIterator>
1623 const InputVector &values,
1624 ForwardIterator local_values_begin,
1625 ForwardIterator local_values_end)
const;
1651 template <
class OutputVector,
typename number>
1654 OutputVector &values)
const;
1687 template <
typename Number>
1755 template <
class OutputVector,
typename number>
1759 OutputVector &values,
1761 const bool perform_check =
false)
const;
1772 template <
class OutputVector,
typename number>
1776 OutputVector &values,
1793 template <
typename number,
typename OutputVector>
1796 OutputVector &global_destination)
const;
1812 template <
typename ForwardIterator,
typename OutputVector>
1815 ForwardIterator local_source_end,
1816 OutputVector &global_destination)
const;
1831 template <
typename ForwardIterator,
typename OutputVector>
1835 ForwardIterator local_source_begin,
1836 ForwardIterator local_source_end,
1837 OutputVector &global_destination)
const;
1845 template <
typename number,
typename OutputMatrix>
1848 OutputMatrix &global_destination)
const;
1854 template <
typename number,
typename OutputMatrix,
typename OutputVector>
1858 OutputMatrix &global_matrix,
1859 OutputVector &global_vector)
const;
1887 std::vector<types::global_dof_index> &dof_indices)
const;
2114 friend struct ::internal::DoFCellAccessorImplementation::Implementation;
2118template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2122 return level_dof_access;
2127template <
int structdim,
int dim,
int spacedim>
2128template <
typename OtherAccessor>
2130 const OtherAccessor &)
2133 ExcMessage(
"You are attempting an illegal conversion between "
2134 "iterator/accessor types. The constructor you call "
2135 "only exists to make certain template constructs "
2136 "easier to write as dimension independent code but "
2137 "the conversion is not valid in the current context."));
2145#include <deal.II/dofs/dof_accessor.templates.h>
DoFAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
typename::internal::DoFHandlerImplementation::Iterators< dim, spacedim, level_dof_access >::line_iterator line(const unsigned int i) const
friend struct ::internal::DoFHandlerImplementation::Policy::Implementation
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
void get_mg_dof_indices(const int level, std::vector< types::global_dof_index > &dof_indices, const types::fe_index fe_index=numbers::invalid_fe_index) const
DoFAccessor< 0, 1, spacedim, level_dof_access > & operator=(const DoFAccessor< 0, 1, spacedim, level_dof_access > &da)=delete
void get_dof_indices(std::vector< types::global_dof_index > &dof_indices, const types::fe_index fe_index=numbers::invalid_fe_index) const
typename::internal::DoFHandlerImplementation::Iterators< 1, spacedim, level_dof_access >::line_iterator line(const unsigned int i) const
DoFAccessor< 0, 1, spacedim, level_dof_access > & operator=(DoFAccessor< 0, 1, spacedim, level_dof_access > &&) noexcept=default
types::global_dof_index vertex_dof_index(const unsigned int vertex, const unsigned int i, const types::fe_index fe_index=numbers::invalid_fe_index) const
static constexpr unsigned int space_dimension
TriaAccessor< 0, 1, spacedim > BaseClass
DoFAccessor< structdim, dim, spacedim, level_dof_access > & operator=(const DoFAccessor< structdim, dim, spacedim, level_dof_access > &da)=delete
types::fe_index nth_active_fe_index(const unsigned int n) const
DoFAccessor(DoFAccessor< structdim, dim, spacedim, level_dof_access > &&)
DoFAccessor(const DoFAccessor< structdim, dim, spacedim, level_dof_access2 > &)
static constexpr unsigned int dimension
TriaIterator< DoFAccessor< structdim, dim, spacedim, level_dof_access > > child(const unsigned int c) const
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
void set_mg_vertex_dof_index(const int level, const unsigned int vertex, const unsigned int i, const types::global_dof_index index, const types::fe_index fe_index=numbers::invalid_fe_index) const
const DoFHandler< 1, spacedim > & get_dof_handler() const
DoFHandler< 1, spacedim > * dof_handler
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 types::fe_index fe_index=numbers::invalid_fe_index) const
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)
typename::internal::DoFHandlerImplementation::Iterators< 1, spacedim, level_dof_access >::quad_iterator quad(const unsigned int i) const
DoFAccessor< structdim, dim, spacedim, level_dof_access > & operator=(DoFAccessor< structdim, dim, spacedim, level_dof_access > &&)
std::set< types::fe_index > get_active_fe_indices() const
DoFAccessor(const DoFAccessor< structdim, dim, spacedim, level_dof_access > &)=default
DoFAccessor(const DoFAccessor< 0, 1, spacedim, level_dof_access > &)=default
void set_mg_dof_indices(const int level, const std::vector< types::global_dof_index > &dof_indices, const types::fe_index fe_index=numbers::invalid_fe_index)
void get_mg_dof_indices(const int level, std::vector< types::global_dof_index > &dof_indices, const types::fe_index fe_index=numbers::invalid_fe_index) const
types::global_dof_index mg_vertex_dof_index(const int level, const unsigned int vertex, const unsigned int i, const types::fe_index fe_index=numbers::invalid_fe_index) const
typename ::internal::DoFAccessorImplementation:: Inheritance< structdim, dimension, space_dimension >::BaseClass BaseClass
DoFAccessor(const Triangulation< dim, spacedim > *tria, const int level, const int index, const DoFHandler< dim, spacedim > *dof_handler)
types::fe_index nth_active_fe_index(const unsigned int n) const
const FiniteElement< 1, spacedim > & get_fe(const types::fe_index fe_index) const
bool operator!=(const DoFAccessor< structdim2, dim2, spacedim2, level_dof_access2 > &) const
void copy_from(const TriaAccessorBase< structdim, dim, spacedim > &da)
types::global_dof_index dof_index(const unsigned int i, const types::fe_index fe_index=numbers::invalid_fe_index) const
typename::internal::DoFHandlerImplementation::Iterators< dim, spacedim, level_dof_access >::quad_iterator quad(const unsigned int i) const
friend class TriaRawIterator
DoFHandler< dimension, space_dimension > AccessorData
void set_dof_index(const unsigned int i, const types::global_dof_index index, const types::fe_index fe_index=numbers::invalid_fe_index) const
void copy_from(const DoFAccessor< structdim, dim, spacedim, level_dof_access2 > &a)
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index fe_index) const
static constexpr unsigned int space_dimension
void get_dof_indices(std::vector< types::global_dof_index > &dof_indices, const types::fe_index fe_index=numbers::invalid_fe_index) const
void set_dof_handler(DoFHandler< 1, spacedim > *dh)
unsigned int n_active_fe_indices() 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
void set_dof_index(const unsigned int i, const types::global_dof_index index, const types::fe_index fe_index=numbers::invalid_fe_index) const
DoFHandler< 1, spacedim > AccessorData
unsigned int n_active_fe_indices() const
TriaIterator< DoFAccessor< 0, 1, spacedim, level_dof_access > > child(const unsigned int c) const
DoFAccessor(const DoFAccessor< structdim2, dim2, spacedim2, level_dof_access2 > &)
static bool is_level_cell()
types::global_dof_index vertex_dof_index(const unsigned int vertex, const unsigned int i, const types::fe_index fe_index=numbers::invalid_fe_index) const
bool fe_index_is_active(const types::fe_index fe_index) const
void copy_from(const DoFAccessor< 0, 1, spacedim, level_dof_access2 > &a)
static constexpr unsigned int dimension
DoFHandler< dim, spacedim > * dof_handler
DoFAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
bool fe_index_is_active(const types::fe_index fe_index) const
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > periodic_neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) 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
DoFHandler< dimension_, space_dimension_ > AccessorData
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
types::fe_index future_fe_index() 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
DoFAccessor< dimension_, dimension_, space_dimension_, level_dof_access > BaseClass
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > child(const unsigned int i) const
TriaIterator< DoFAccessor< dimension_ - 1, dimension_, space_dimension_, level_dof_access > > face_iterator
DoFCellAccessor(const Triangulation< dimension_, space_dimension_ > *tria, const int level, const int index, const AccessorData *local_data)
void set_future_fe_index(const types::fe_index i) const
void distribute_local_to_global(const Vector< number > &local_source, OutputVector &global_destination) const
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > neighbor_or_periodic_neighbor(const unsigned int i) 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
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
static const unsigned int dim
void set_active_fe_index(const types::fe_index i) const
void distribute_local_to_global_by_interpolation(const Vector< number > &local_values, OutputVector &values, const types::fe_index fe_index=numbers::invalid_fe_index) const
DoFCellAccessor(const DoFAccessor< structdim2, dim2, spacedim2, level_dof_access2 > &)
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > periodic_neighbor(const unsigned int i) const
DoFCellAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
static const unsigned int spacedim
void get_interpolated_dof_values(const ReadVector< Number > &values, Vector< Number > &interpolated_values, const types::fe_index fe_index=numbers::invalid_fe_index) const
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 set_dof_indices(const std::vector< types::global_dof_index > &dof_indices)
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_dof_values(const Vector< number > &local_values, OutputVector &values) const
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
void set_dof_values_by_interpolation(const Vector< number > &local_values, OutputVector &values, const types::fe_index fe_index=numbers::invalid_fe_index, const bool perform_check=false) 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
void get_dof_values(const ReadVector< Number > &values, ForwardIterator local_values_begin, ForwardIterator local_values_end) const
~DoFCellAccessor()=default
DoFHandler< dimension_, space_dimension_ > Container
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
types::fe_index active_fe_index() const
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > parent() const
static const types::fe_index default_fe_index
typename InvalidAccessor< structdim, dim, spacedim >::AccessorData AccessorData
DoFInvalidAccessor(const void *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
InvalidAccessor(const void *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
friend class TriaIterator
TriaAccessorBase(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *=nullptr)
TriaAccessor(const Triangulation< 1, spacedim > *tria, const VertexKind vertex_kind, const unsigned int vertex_index)
unsigned int vertex_index(const unsigned int i=0) const
Point< spacedim > & vertex(const unsigned int i=0) const
const Triangulation< 1, spacedim > * tria
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcVectorNotEmpty()
#define DeclException0(Exception0)
static ::ExceptionBase & ExcVectorDoesNotMatch()
static ::ExceptionBase & ExcInvalidObject()
static ::ExceptionBase & ExcVectorNotEmpty()
static ::ExceptionBase & ExcMatrixDoesNotMatch()
#define Assert(cond, exc)
static ::ExceptionBase & ExcCantCompareIterators()
static ::ExceptionBase & ExcNotActive()
static ::ExceptionBase & ExcMatrixDoesNotMatch()
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcInvalidObject()
static ::ExceptionBase & ExcNotActive()
static ::ExceptionBase & ExcCantCompareIterators()
static ::ExceptionBase & ExcVectorDoesNotMatch()
static ::ExceptionBase & ExcMessage(std::string arg1)
const types::fe_index invalid_fe_index
unsigned short int fe_index
unsigned int global_dof_index
static constexpr unsigned int faces_per_cell
static constexpr unsigned int max_children_per_cell
::TriaAccessor< structdim, dim, spacedim > BaseClass
::CellAccessor< dim, spacedim > BaseClass