![]() |
Reference documentation for deal.II version 8.1.0
|
#include <dof_accessor.h>
Public Types | |
typedef typename::internal::DoFAccessor::Inheritance< structdim, dimension, space_dimension >::BaseClass | BaseClass |
typedef DH | AccessorData |
![]() | |
typedef TriaAccessorBase< structdim, dim, spacedim >::AccessorData | AccessorData |
![]() | |
typedef void * | LocalData |
Public Member Functions | |
const DH & | get_dof_handler () const |
template<bool level_dof_access2> | |
void | copy_from (const DoFAccessor< structdim, DH, level_dof_access2 > &a) |
void | copy_from (const TriaAccessorBase< structdim, DH::dimension, DH::space_dimension > &da) |
TriaIterator< DoFAccessor< structdim, DH, level_dof_access > > | parent () const |
DeclException0 (ExcInvalidObject) | |
DeclException0 (ExcVectorNotEmpty) | |
DeclException0 (ExcVectorDoesNotMatch) | |
DeclException0 (ExcMatrixDoesNotMatch) | |
DeclException0 (ExcNotActive) | |
DeclException0 (ExcCantCompareIterators) | |
template<int dim2, class DH2 , bool level_dof_access2> | |
bool | operator== (const DoFAccessor< dim2, DH2, level_dof_access2 > &) const |
template<int dim2, class DH2 , bool level_dof_access2> | |
bool | operator!= (const DoFAccessor< dim2, DH2, level_dof_access2 > &) const |
template<template< int, int > class DH, int spacedim, bool lda> | |
DoFAccessor (const Triangulation< 1, spacedim > *tria, const typename TriaAccessor< 0, 1, spacedim >::VertexKind vertex_kind, const unsigned int vertex_index, const DH< 1, spacedim > *dof_handler) | |
template<template< int, int > class DH, int spacedim, bool lda> | |
DoFAccessor (const Triangulation< 1, spacedim > *, const int, const int, const DH< 1, spacedim > *) | |
template<bool lda2> | |
void | copy_from (const DoFAccessor< 0, DH< 1, spacedim >, lda2 > &a) |
Constructors | |
DoFAccessor () | |
DoFAccessor (const Triangulation< DH::dimension, DH::space_dimension > *tria, const int level, const int index, const DH *local_data) | |
template<int structdim2, int dim2, int spacedim2> | |
DoFAccessor (const InvalidAccessor< structdim2, dim2, spacedim2 > &) | |
template<int dim2, class DH2 , bool level_dof_access2> | |
DoFAccessor (const DoFAccessor< dim2, DH2, level_dof_access2 > &) | |
template<bool level_dof_access2> | |
DoFAccessor (const DoFAccessor< structdim, DH, level_dof_access2 > &) | |
Accessing sub-objects | |
TriaIterator< DoFAccessor< structdim, DH, level_dof_access > > | child (const unsigned int c) const |
typename::internal::DoFHandler::Iterators< DH, level_dof_access >::line_iterator | line (const unsigned int i) const |
typename::internal::DoFHandler::Iterators< DH, level_dof_access >::quad_iterator | quad (const unsigned int i) const |
Accessing the DoF indices of this object | |
void | get_dof_indices (std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DH::default_fe_index) const |
void | get_mg_dof_indices (const int level, std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DH::default_fe_index) const |
void | set_mg_dof_indices (const int level, const std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DH::default_fe_index) |
types::global_dof_index | vertex_dof_index (const unsigned int vertex, const unsigned int i, const unsigned int fe_index=DH::default_fe_index) 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=DH::default_fe_index) const |
types::global_dof_index | dof_index (const unsigned int i, const unsigned int fe_index=DH::default_fe_index) const |
types::global_dof_index | mg_dof_index (const int level, const unsigned int i) const |
Accessing the finite element associated with this object | |
unsigned int | n_active_fe_indices () const |
unsigned int | nth_active_fe_index (const unsigned int n) const |
bool | fe_index_is_active (const unsigned int fe_index) const |
const FiniteElement< DH::dimension, DH::space_dimension > & | get_fe (const unsigned int fe_index) const |
![]() | |
TriaAccessor (const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *local_data=0) | |
template<int structdim2, int dim2, int spacedim2> | |
TriaAccessor (const InvalidAccessor< structdim2, dim2, spacedim2 > &) | |
template<int structdim2, int dim2, int spacedim2> | |
TriaAccessor (const TriaAccessor< structdim2, dim2, spacedim2 > &) | |
bool | used () const |
int | parent_index () const |
unsigned int | vertex_index (const unsigned int i) const |
Point< spacedim > & | vertex (const unsigned int i) const |
typename::internal::Triangulation::Iterators< dim, spacedim >::line_iterator | line (const unsigned int i) const |
unsigned int | line_index (const unsigned int i) const |
typename::internal::Triangulation::Iterators< dim, spacedim >::quad_iterator | quad (const unsigned int i) const |
unsigned int | quad_index (const unsigned int i) const |
bool | face_orientation (const unsigned int face) const |
bool | face_flip (const unsigned int face) const |
bool | face_rotation (const unsigned int face) const |
bool | line_orientation (const unsigned int line) const |
bool | has_children () const |
unsigned int | n_children () const |
unsigned int | number_of_children () const |
unsigned int | max_refinement_depth () const |
TriaIterator< TriaAccessor< structdim, dim, spacedim > > | child (const unsigned int i) const |
TriaIterator< TriaAccessor< structdim, dim, spacedim > > | isotropic_child (const unsigned int i) const |
RefinementCase< structdim > | refinement_case () const |
int | child_index (const unsigned int i) const |
int | isotropic_child_index (const unsigned int i) const |
types::boundary_id | boundary_indicator () const |
void | set_boundary_indicator (const types::boundary_id) const |
void | set_all_boundary_indicators (const types::boundary_id) const |
bool | at_boundary () const |
const Boundary< dim, spacedim > & | get_boundary () const |
bool | user_flag_set () const |
void | set_user_flag () const |
void | clear_user_flag () const |
void | recursively_set_user_flag () const |
void | recursively_clear_user_flag () const |
void | clear_user_data () const |
void | set_user_pointer (void *p) const |
void | clear_user_pointer () const |
void * | user_pointer () const |
void | recursively_set_user_pointer (void *p) const |
void | recursively_clear_user_pointer () const |
void | set_user_index (const unsigned int p) const |
void | clear_user_index () const |
unsigned int | user_index () const |
void | recursively_set_user_index (const unsigned int p) const |
void | recursively_clear_user_index () const |
double | diameter () const |
double | extent_in_direction (const unsigned int axis) const |
double | minimum_vertex_distance () const |
Point< spacedim > | center () const |
Point< spacedim > | barycenter () const |
double | measure () const |
bool | is_translation_of (const TriaIterator< TriaAccessor< structdim, dim, spacedim > > &o) const |
![]() | |
int | level () const |
int | index () const |
IteratorState::IteratorStates | state () const |
const Triangulation< dim, spacedim > & | get_triangulation () const |
Static Public Member Functions | |
static bool | is_level_cell () |
Static Public Attributes | |
static const unsigned int | dimension =DH::dimension |
static const unsigned int | space_dimension =DH::space_dimension |
![]() | |
static const unsigned int | space_dimension = spacedim |
static const unsigned int | dimension = dim |
static const unsigned int | structure_dimension = structdim |
Protected Member Functions | |
void | set_dof_handler (DH *dh) |
void | set_dof_index (const unsigned int i, const types::global_dof_index index, const unsigned int fe_index=DH::default_fe_index) const |
void | set_mg_dof_index (const int level, const unsigned int i, const types::global_dof_index index) const |
void | set_vertex_dof_index (const unsigned int vertex, const unsigned int i, const types::global_dof_index index, const unsigned int fe_index=DH::default_fe_index) 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=DH::default_fe_index) const |
![]() | |
TriaAccessorBase (const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *=0) | |
TriaAccessorBase (const TriaAccessorBase &) | |
void | copy_from (const TriaAccessorBase &) |
TriaAccessorBase & | operator= (const TriaAccessorBase &) |
bool | operator< (const TriaAccessorBase &other) const |
void | operator= (const TriaAccessorBase *) |
bool | operator== (const TriaAccessorBase &) const |
bool | operator!= (const TriaAccessorBase &) const |
::internal::Triangulation::TriaObjects<::internal::Triangulation::TriaObject< structdim > > & | objects () const |
void | operator++ () |
void | operator-- () |
Protected Attributes | |
DH * | dof_handler |
![]() | |
typename::internal::TriaAccessor::PresentLevelType< structdim, dim >::type | present_level |
int | present_index |
const Triangulation< dim, spacedim > * | tria |
Private Member Functions | |
DoFAccessor< structdim, DH, level_dof_access > & | operator= (const DoFAccessor< structdim, DH, level_dof_access > &da) |
Friends | |
template<typename > | |
class | TriaRawIterator |
template<int , class , bool > | |
class | DoFAccessor |
template<int dim, int spacedim> | |
class | DoFHandler |
template<int dim, int spacedim> | |
class | hp::DoFHandler |
struct | ::internal::DoFHandler::Policy::Implementation |
struct | ::internal::DoFHandler::Implementation |
struct | ::internal::hp::DoFHandler::Implementation |
struct | ::internal::DoFCellAccessor::Implementation |
struct | ::internal::DoFAccessor::Implementation |
Additional Inherited Members | |
![]() | |
typedef void | AccessorData |
A class that gives access to the degrees of freedom stored in a DoFHandler, MGSoDHandler, or hp::DoFHandler object. Accessors are used to, access the data that pertains to edges, faces, and cells of a triangulation. The concept is explained in more detail in connection to Iterators on mesh-like containers.
This class follows mainly the route laid out by the accessor library declared in the triangulation library (TriaAccessor). It enables the user to access the degrees of freedom on lines, quads, or hexes. The first template argument of this class determines the dimensionality of the object under consideration: 1 for lines, 2 for quads, and 3 for hexes. The second argument denotes the type of DoF handler we should work on. It can either be DoFHandler or hp::DoFHandler. From the second template argument we also deduce the dimension of the Triangulation this object refers to as well as the dimension of the space into which it is embedded. Finally, the template argument level_dof_access
governs the behavior of the function get_active_or_mg_dof_indices(). See the section on Generic loops below.
Usage is best to happen through the typedefs to the various kinds of iterators provided by the DoFHandler and hp::DoFHandler classes, since they are more secure to changes in the class naming and template interface as well as providing easier typing (much less complicated names!).
Many loops look very similar, whether they operate on the active dofs of the active cells of the Triangulation or on the level dodfs of a single level or the whole grid hierarchy. In order to use polymorphism in such loops, they access degrees of freedom through the function get_active_or_mg_dof_indices(), which changes behavior according to the third template argument. If the argument is false, then the active dofs of active cells are accessed. If it is true, the level dofs are used. DoFHandler has functions, for instance begin() and begin_mg(), which return either type or the other. Additionally, they can be cast into each other, in case this is needed, since they access the same data.
It is highly recommended to use the function get_active_or_mg_dof_indices() in generic loops in lieu of get_dof_indices() or get_mg_dof_indices().
If the structural dimension given by the first template argument equals the dimension of the DoFHandler (given as the second template argument), then we are obviously dealing with cells, rather than lower-dimensional objects. In that case, inheritance is from CellAccessor, to provide access to all the cell specific information afforded by that class. Otherwise, i.e. for lower-dimensional objects, inheritance is from TriaAccessor.
There is a DoFCellAccessor class that provides the equivalent to the CellAccessor class.
Definition at line 189 of file dof_accessor.h.
typedef typename ::internal::DoFAccessor::Inheritance<structdim, dimension, space_dimension>::BaseClass DoFAccessor< structdim, DH, level_dof_access >::BaseClass |
Declare a typedef to the base class to make accessing some of the exception classes simpler.
Definition at line 215 of file dof_accessor.h.
typedef DH DoFAccessor< structdim, DH, level_dof_access >::AccessorData |
Data type passed by the iterator class.
Definition at line 220 of file dof_accessor.h.
|
inline |
Default constructor. Provides an accessor that can't be used.
Definition at line 43 of file dof_accessor.templates.h.
|
inline |
Constructor
Definition at line 52 of file dof_accessor.templates.h.
DoFAccessor< DH, spacedim, lda >::DoFAccessor | ( | const InvalidAccessor< structdim2, dim2, spacedim2 > & | ) |
Conversion constructor. This constructor exists to make certain constructs simpler to write in dimension independent code. For example, it allows assigning a face iterator to a line iterator, an operation that is useful in 2d but doesn't make any sense in 3d. The constructor here exists for the purpose of making the code conform to C++ but it will unconditionally abort; in other words, assigning a face iterator to a line iterator is better put into an if-statement that checks that the dimension is two, and assign to a quad iterator in 3d (an operator that, without this constructor would be illegal if we happen to compile for 2d).
Definition at line 69 of file dof_accessor.templates.h.
|
inline |
Another conversion operator between objects that don't make sense, just like the previous one.
Definition at line 79 of file dof_accessor.templates.h.
|
inline |
Copy constructor allowing to switch level access and active access.
Definition at line 90 of file dof_accessor.templates.h.
|
inline |
Return a handle on the DoFHandler object which we are using.
Definition at line 111 of file dof_accessor.templates.h.
|
inline |
Implement the copy operator needed for the iterator classes.
Definition at line 134 of file dof_accessor.templates.h.
|
inline |
Copy operator used by the iterator class. Keeps the previously set dof handler, but sets the object coordinates of the TriaAccessor.
Definition at line 121 of file dof_accessor.templates.h.
|
inlinestatic |
Tell the caller whether get_active_or_mg_dof_indices() accesses active or level dofs.
Definition at line 2205 of file dof_accessor.h.
|
inline |
Return an iterator pointing to the the parent.
Definition at line 187 of file dof_accessor.templates.h.
|
inline |
Return an iterator pointing to the the c-th
child.
Definition at line 171 of file dof_accessor.templates.h.
|
inline |
Pointer to the ith
line bounding this object. If the current object is a line itself, then the only valid index is i
equals to zero, and the function returns an iterator to itself.
Definition at line 2117 of file dof_accessor.templates.h.
|
inline |
Pointer to the ith
quad bounding this object. If the current object is a quad itself, then the only valid index is i
equals to zero, and the function returns an iterator to itself.
Definition at line 2150 of file dof_accessor.templates.h.
|
inline |
Return the global indices of the degrees of freedom located on this object in the standard ordering defined by the finite element (i.e., dofs on vertex 0, dofs on vertex 1, etc, dofs on line 0, dofs on line 1, etc, dofs on quad 0, etc.) This function is only available on active objects (see this glossary entry).
The cells needs to be an active cell (and not artificial in a parallel distributed computation).
The vector has to have the right size before being passed to this function.
The last argument denotes the finite element index. For the standard DoFHandler class, this value must be equal to its default value since that class only supports the same finite element on all cells anyway.
However, for hp objects (i.e. the hp::DoFHandler class), different finite element objects may be used on different cells. On faces between two cells, as well as vertices, there may therefore be two sets of degrees of freedom, one for each of the finite elements used on the adjacent cells. In order to specify which set of degrees of freedom to work on, the last argument is used to disambiguate. Finally, if this function is called for a cell object, there can only be a single set of degrees of freedom, and fe_index has to match the result of active_fe_index().
For cells, there is only a single possible finite element index (namely the one for that cell, returned by cell->active_fe_index
. Consequently, the derived DoFCellAccessor class has an overloaded version of this function that calls the present function with cell->active_fe_index
as last argument.
Definition at line 1889 of file dof_accessor.templates.h.
|
inline |
Return the global multilevel indices of the degrees of freedom that live on the current object with respect to the given level within the multigrid hierarchy. The indices refer to the local numbering for the level this line lives on.
Definition at line 1951 of file dof_accessor.templates.h.
|
inline |
Sets the level DoF indices that are returned by get_mg_dof_indices.
Definition at line 2003 of file dof_accessor.templates.h.
|
inline |
Global DoF index of the i degree associated with the vertexth
vertex of the present cell.
The last argument denotes the finite element index. For the standard DoFHandler class, this value must be equal to its default value since that class only supports the same finite element on all cells anyway.
However, for hp objects (i.e. the hp::DoFHandler class), different finite element objects may be used on different cells. On faces between two cells, as well as vertices, there may therefore be two sets of degrees of freedom, one for each of the finite elements used on the adjacent cells. In order to specify which set of degrees of freedom to work on, the last argument is used to disambiguate. Finally, if this function is called for a cell object, there can only be a single set of degrees of freedom, and fe_index has to match the result of active_fe_index().
Definition at line 1569 of file dof_accessor.templates.h.
|
inline |
Returns the global DoF index of the i
th degree of freedom associated with the vertex
th vertex on level level
. Also see vertex_dof_index().
Definition at line 1585 of file dof_accessor.templates.h.
|
inline |
Index of the ith degree of freedom of this object.
The last argument denotes the finite element index. For the standard DoFHandler class, this value must be equal to its default value since that class only supports the same finite element on all cells anyway.
However, for hp objects (i.e. the hp::DoFHandler class), different finite element objects may be used on different cells. On faces between two cells, as well as vertices, there may therefore be two sets of degrees of freedom, one for each of the finite elements used on the adjacent cells. In order to specify which set of degrees of freedom to work on, the last argument is used to disambiguate. Finally, if this function is called for a cell object, there can only be a single set of degrees of freedom, and fe_index has to match the result of active_fe_index().
Definition at line 1474 of file dof_accessor.templates.h.
|
inline |
Returns the dof_index on the given level. Also see dof_index.
Definition at line 1490 of file dof_accessor.templates.h.
|
inline |
Return the number of finite elements that are active on a given object.
For non-hp DoFHandler objects, the answer is of course always one. However, for hp::DoFHandler objects, this isn't the case: If this is a cell, the answer is of course one. If it is a face, the answer may be one or two, depending on whether the two adjacent cells use the same finite element or not. If it is an edge in 3d, the possible return value may be one or any other value larger than that.
Definition at line 1519 of file dof_accessor.templates.h.
|
inline |
Return the n-th
active fe index on this object. For cells and all non-hp objects, there is only a single active fe index, so the argument must be equal to zero. For lower-dimensional hp objects, there are n_active_fe_indices() active finite elements, and this function can be queried for their indices.
Definition at line 1535 of file dof_accessor.templates.h.
|
inline |
Return true if the finite element with given index is active on the present object. For non-hp DoF accessors, this is of course the case only if fe_index
equals zero. For cells, it is the case if fe_index
equals active_fe_index() of this cell. For faces and other lower-dimensional objects, there may be more than one fe_index
that are active on any given object (see n_active_fe_indices()).
Definition at line 1552 of file dof_accessor.templates.h.
|
inline |
Return a reference to the finite element used on this object with the given fe_index
. fe_index
must be used on this object, i.e. fe_index_is_active(fe_index)
must return true.
Definition at line 1673 of file dof_accessor.templates.h.
|
inline |
Compare for equality. Return true
if the two accessors refer to the same object.
The template parameters of this function allow for a comparison of very different objects. Therefore, some of them are disabled. Namely, if the dimension, or the dof handler of the two objects differ, an exception is generated. It can be expected that this is an unwanted comparison.
The template parameter level_dof_access2
is ignored, such that an iterator with level access can be equal to one with access to the active degrees of freedom.
Definition at line 146 of file dof_accessor.templates.h.
|
inline |
Compare for inequality. The boolean not of operator==().
Definition at line 159 of file dof_accessor.templates.h.
|
inlineprotected |
Reset the DoF handler pointer.
Definition at line 100 of file dof_accessor.templates.h.
|
inlineprotected |
Set the index of the ith degree of freedom of this object to index
.
The last argument denotes the finite element index. For the standard DoFHandler class, this value must be equal to its default value since that class only supports the same finite element on all cells anyway.
However, for hp objects (i.e. the hp::DoFHandler class), different finite element objects may be used on different cells. On faces between two cells, as well as vertices, there may therefore be two sets of degrees of freedom, one for each of the finite elements used on the adjacent cells. In order to specify which set of degrees of freedom to work on, the last argument is used to disambiguate. Finally, if this function is called for a cell object, there can only be a single set of degrees of freedom, and fe_index has to match the result of active_fe_index().
Definition at line 1500 of file dof_accessor.templates.h.
|
inlineprotected |
Set the global index of the i degree on the vertex-th
vertex of the present cell to index
.
The last argument denotes the finite element index. For the standard DoFHandler class, this value must be equal to its default value since that class only supports the same finite element on all cells anyway.
However, for hp objects (i.e. the hp::DoFHandler class), different finite element objects may be used on different cells. On faces between two cells, as well as vertices, there may therefore be two sets of degrees of freedom, one for each of the finite elements used on the adjacent cells. In order to specify which set of degrees of freedom to work on, the last argument is used to disambiguate. Finally, if this function is called for a cell object, there can only be a single set of degrees of freedom, and fe_index has to match the result of active_fe_index().
Definition at line 1601 of file dof_accessor.templates.h.
|
private |
Copy operator. This is normally used in a context like iterator a,b; *a=*b;
. Presumably, the intent here is to copy the object pointed to by b
to the object pointed to by a
. However, the result of dereferencing an iterator is not an object but an accessor; consequently, this operation is not useful for iterators on triangulations. We declare this function here private, thus it may not be used from outside. Furthermore it is not implemented and will give a linker error if used anyway.
|
friend |
Iterator classes need to be friends because they need to access operator== and operator!=.
Definition at line 806 of file dof_accessor.h.
|
friend |
Make the DoFHandler class a friend so that it can call the set_xxx() functions.
Definition at line 835 of file dof_accessor.h.
|
static |
A static variable that allows users of this class to discover the value of the second template argument.
Definition at line 198 of file dof_accessor.h.
|
static |
A static variable that allows users of this class to discover the value of the third template argument.
Definition at line 205 of file dof_accessor.h.
|
protected |
Store the address of the DoFHandler object to be accessed.
Definition at line 675 of file dof_accessor.h.