Classes | |
class | ExcCantCompareIterators |
class | ExcInvalidObject |
class | ExcMatrixDoesNotMatch |
class | ExcNotActive |
class | ExcVectorDoesNotMatch |
class | ExcVectorNotEmpty |
Public Types | |
typedef internal::DoFAccessor::Inheritance < structdim, dimension, space_dimension >::BaseClass | BaseClass |
typedef DH | AccessorData |
Public Member Functions | |
const DH & | get_dof_handler () const |
DoFAccessor< structdim, DH > & | operator= (const DoFAccessor< structdim, DH > &da) |
void | copy_from (const DoFAccessor< structdim, DH > &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 > | |
DoFAccessor (const DoFAccessor< dim2, DH2 > &) | |
Accessing sub-objects | |
TriaIterator< DoFAccessor < structdim, DH > > | child (const unsigned int c) const |
internal::DoFHandler::Iterators < DH >::line_iterator | line (const unsigned int i) const |
internal::DoFHandler::Iterators < DH >::quad_iterator | quad (const unsigned int i) const |
Accessing the DoF indices of this object | |
void | get_dof_indices (std::vector< unsigned int > &dof_indices, const unsigned int fe_index=DH::default_fe_index) const |
unsigned int | vertex_dof_index (const unsigned int vertex, const unsigned int i, const unsigned int fe_index=DH::default_fe_index) const |
unsigned int | dof_index (const unsigned int i, const unsigned int fe_index=DH::default_fe_index) 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 |
Static Public Attributes | |
static const unsigned int | dimension = DH::dimension |
static const unsigned int | space_dimension = DH::space_dimension |
Protected Member Functions | |
bool | operator== (const DoFAccessor &) const |
bool | operator!= (const DoFAccessor &) const |
void | set_dof_handler (DH *dh) |
void | set_dof_index (const unsigned int i, const unsigned int index, const unsigned int fe_index=DH::default_fe_index) const |
void | set_vertex_dof_index (const unsigned int vertex, const unsigned int i, const unsigned int index, const unsigned int fe_index=DH::default_fe_index) const |
Protected Attributes | |
DH * | dof_handler |
Friends | |
class | TriaRawIterator |
class | DoFHandler |
class | hp::DoFHandler |
class | internal::DoFHandler::Implementation |
class | internal::hp::DoFHandler::Implementation |
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. From the second template argument we can deduce the dimensionality of the triangulation to which this object belongs as well as the dimensionality of the space in which it is embedded. The second argument also denotes the type of DoF handler we should work on. It can either be DoFHandler or hp::DoFHandler.
Depending on whether the structural dimension of the object accessed equals the dimension on which the DoF handler object operates, this class is derived from CellAccessor or TriaAccessor. This means that, for example accessors to quads in 2d have access to all the mesh aspects of cells, whereas accessors to quads in 3d can only access things that make sense for faces.
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!).
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.
typedef internal::DoFAccessor::Inheritance<structdim, dimension, space_dimension>::BaseClass DoFAccessor< structdim, DH >::BaseClass |
Declare a typedef to the base class to make accessing some of the exception classes simpler.
Reimplemented in DoFCellAccessor< DH >, MGDoFAccessor< structdim, dim, spacedim >, MGDoFCellAccessor< dim, spacedim >, and MGDoFAccessor< dim, dim, spacedim >.
typedef DH DoFAccessor< structdim, DH >::AccessorData |
Data type passed by the iterator class.
Reimplemented from TriaAccessor< celldim, dim, spacedim >.
Reimplemented in DoFCellAccessor< DH >, MGDoFAccessor< structdim, dim, spacedim >, MGDoFCellAccessor< dim, spacedim >, and MGDoFAccessor< dim, dim, spacedim >.
DoFAccessor< structdim, DH >::DoFAccessor | ( | ) |
Default constructor. Provides an accessor that can't be used.
DoFAccessor< structdim, DH >::DoFAccessor | ( | const Triangulation< DH::dimension, DH::space_dimension > * | tria, | |
const int | level, | |||
const int | index, | |||
const DH * | local_data | |||
) |
Constructor
DoFAccessor< structdim, DH >::DoFAccessor | ( | const InvalidAccessor< structdim2, dim2, spacedim2 > & | ) | [inline] |
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).
DoFAccessor< structdim, DH >::DoFAccessor | ( | const DoFAccessor< dim2, DH2 > & | ) | [inline] |
Another conversion operator between objects that don't make sense, just like the previous one.
const DH& DoFAccessor< structdim, DH >::get_dof_handler | ( | ) | const |
Return a handle on the DoFHandler object which we are using.
DoFAccessor<structdim,DH>& DoFAccessor< structdim, DH >::operator= | ( | const DoFAccessor< structdim, DH > & | da | ) |
Copy operator.
void DoFAccessor< structdim, DH >::copy_from | ( | const DoFAccessor< structdim, DH > & | a | ) |
Implement the copy operator needed for the iterator classes.
TriaIterator<DoFAccessor<structdim,DH> > DoFAccessor< structdim, DH >::child | ( | const unsigned int | c | ) | const |
Return an iterator pointing to the the c-th
child.
Reimplemented from TriaAccessor< celldim, dim, spacedim >.
Reimplemented in DoFCellAccessor< DH >, MGDoFAccessor< structdim, dim, spacedim >, MGDoFCellAccessor< dim, spacedim >, and MGDoFAccessor< dim, dim, spacedim >.
internal::DoFHandler::Iterators<DH>::line_iterator DoFAccessor< structdim, DH >::line | ( | const unsigned int | i | ) | const |
Pointer to the ith
line bounding this object.
Reimplemented from TriaAccessor< celldim, dim, spacedim >.
Reimplemented in MGDoFAccessor< structdim, dim, spacedim >, and MGDoFAccessor< dim, dim, spacedim >.
internal::DoFHandler::Iterators<DH>::quad_iterator DoFAccessor< structdim, DH >::quad | ( | const unsigned int | i | ) | const |
Pointer to the ith
quad bounding this object.
Reimplemented from TriaAccessor< celldim, dim, spacedim >.
Reimplemented in MGDoFAccessor< structdim, dim, spacedim >, and MGDoFAccessor< dim, dim, spacedim >.
void DoFAccessor< structdim, DH >::get_dof_indices | ( | std::vector< unsigned int > & | dof_indices, | |
const unsigned int | fe_index = DH::default_fe_index | |||
) | const |
Return the indices of the dofs of this object in the standard ordering: dofs on vertex 0, dofs on vertex 1, etc, dofs on line 0, dofs on line 1, etc, dofs on quad 0, etc.
The vector has to have the right size before being passed to this function.
This function is most often used on active objects (edges, faces, cells). It can be used on non-active objects as well (i.e. objects that have children), but only if the finite element under consideration has degrees of freedom exclusively on vertices. Otherwise, the function doesn't make much sense, since for example inactive edges do not have degrees of freedom associated with them at all.
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.
unsigned int DoFAccessor< structdim, DH >::vertex_dof_index | ( | const unsigned int | vertex, | |
const unsigned int | i, | |||
const unsigned int | fe_index = DH::default_fe_index | |||
) | const |
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().
unsigned int DoFAccessor< structdim, DH >::dof_index | ( | const unsigned int | i, | |
const unsigned int | fe_index = DH::default_fe_index | |||
) | const |
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().
unsigned int DoFAccessor< structdim, DH >::n_active_fe_indices | ( | ) | const |
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.
unsigned int DoFAccessor< structdim, DH >::nth_active_fe_index | ( | const unsigned int | n | ) | const |
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.
bool DoFAccessor< structdim, DH >::fe_index_is_active | ( | const unsigned int | fe_index | ) | const |
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()).
const FiniteElement<DH::dimension,DH::space_dimension>& DoFAccessor< structdim, DH >::get_fe | ( | const unsigned int | fe_index | ) | const |
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.
bool DoFAccessor< structdim, DH >::operator== | ( | const DoFAccessor< structdim, DH > & | ) | const [protected] |
Compare for equality.
bool DoFAccessor< structdim, DH >::operator!= | ( | const DoFAccessor< structdim, DH > & | ) | const [protected] |
Compare for inequality.
void DoFAccessor< structdim, DH >::set_dof_handler | ( | DH * | dh | ) | [protected] |
Reset the DoF handler pointer.
void DoFAccessor< structdim, DH >::set_dof_index | ( | const unsigned int | i, | |
const unsigned int | index, | |||
const unsigned int | fe_index = DH::default_fe_index | |||
) | const [protected] |
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().
void DoFAccessor< structdim, DH >::set_vertex_dof_index | ( | const unsigned int | vertex, | |
const unsigned int | i, | |||
const unsigned int | index, | |||
const unsigned int | fe_index = DH::default_fe_index | |||
) | const [protected] |
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().
friend class TriaRawIterator [friend] |
Iterator classes need to be friends because they need to access operator== and operator!=.
Reimplemented from TriaAccessorBase< structdim, dim, spacedim >.
friend class DoFHandler [friend] |
Make the DoFHandler class a friend so that it can call the set_xxx() functions.
Reimplemented in DoFCellAccessor< DH >.
friend class hp::DoFHandler [friend] |
friend class internal::DoFHandler::Implementation [friend] |
friend class internal::hp::DoFHandler::Implementation [friend] |
const unsigned int DoFAccessor< structdim, DH >::dimension = DH::dimension [static] |
A static variable that allows users of this class to discover the value of the second template argument.
Reimplemented from TriaAccessorBase< structdim, dim, spacedim >.
const unsigned int DoFAccessor< structdim, DH >::space_dimension = DH::space_dimension [static] |
A static variable that allows users of this class to discover the value of the third template argument.
Reimplemented from TriaAccessorBase< structdim, dim, spacedim >.
DH* DoFAccessor< structdim, DH >::dof_handler [protected] |
Store the address of the DoFHandler object to be accessed.