DoFAccessor< structdim, DH > Class Template Reference
[Degrees of FreedomAccessor classes of the mesh iterators]

Inheritance diagram for DoFAccessor< structdim, DH >:

Inheritance graph
[legend]

List of all members.

Public Types

typedef
DoFObjectAccessor_Inheritance
< structdim, DH::dimension >
::BaseClass 
BaseClass
typedef DH AccessorData

Public Member Functions

 DoFAccessor ()
 DoFAccessor (const Triangulation< DH::dimension > *tria, const int level, const int index, const DH *local_data)
void set_dof_handler (DH *dh)
const DH & get_dof_handler () const
DoFAccessor< structdim, DH > & operator= (const DoFAccessor< structdim, DH > &da)
void copy_from (const DoFAccessor< structdim, DH > &a)
TriaIterator< DH::dimension,
DoFObjectAccessor< structdim,
DH > > 
child (const unsigned int c) const
unsigned int vertex_dof_index (const unsigned int vertex, const unsigned int i, 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
unsigned int dof_index (const unsigned int i, const unsigned int fe_index=DH::default_fe_index) const
void set_dof_index (const unsigned int i, const unsigned int index, const unsigned int fe_index=DH::default_fe_index) const
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 > & 
get_fe (const unsigned int fe_index) const

Protected Member Functions

bool operator== (const DoFAccessor &) const
bool operator!= (const DoFAccessor &) const

Protected Attributes

DH * dof_handler

Friends

class TriaRawIterator

Classes

class  ExcCantCompareIterators
class  ExcInvalidObject
class  ExcMatrixDoesNotMatch
class  ExcNotActive
class  ExcVectorDoesNotMatch
class  ExcVectorNotEmpty


Detailed Description

template<int structdim, class DH>
class DoFAccessor< structdim, DH >

Define the base class for accessors to the degrees of freedom. Accessors are used to, well, 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 the lines (there are similar versions for the DoFs on quads, etc), where the dimension of the underlying triangulation does not really matter (i.e. this accessor works with the lines in 1D-, 2D-, etc dimensions).

The first template argument of this class refers to the structural dimension of the thing accessed: it is 1 for lines, 2 for quads, or 3 for hexes. The second argument denotes the type of DoF handler we should work on. It can either be DoFHandler<dim> or hp::DoFHandler<dim>. The space dimension of the object we are to work on is automatically extracted from the DH template argument.

Depending on whether the structural dimension of the object accessed equals the space dimension on which the DoF handler object operates, this class is derived from CellAccessor or TriaObjectAccessor. 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

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!).

Inheritance

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 TriaObjectAccessor.

There are classes DoFObjectAccessor<1>, DoFObjectAccessor<2>, and DoFObjectAccessor<3> derived from the present class that provide things that are specific to objects of a particular dimensionality. For example, DoFObjectAccessor<2> allows access to the DoF iterators pointing to the lines that bound that quadrilateral, which is obviously something that a line can't do. Consequently, these few operations are declared in these derived classes.

In addition, there is a DoFCellAccessor class that provides the equivalent to the CellAccessor class.

Author:
Wolfgang Bangerth, 1998, 2006

Member Typedef Documentation

template<int structdim, class DH>
typedef DoFObjectAccessor_Inheritance<structdim, DH::dimension>::BaseClass DoFAccessor< structdim, DH >::BaseClass

Declare a typedef to the base class to make accessing some of the exception classes simpler.

Reimplemented in DoFObjectAccessor< 1, DH >, DoFObjectAccessor< 2, DH >, and DoFObjectAccessor< 3, DH >.

template<int structdim, class DH>
typedef DH DoFAccessor< structdim, DH >::AccessorData


Constructor & Destructor Documentation

template<int structdim, class DH>
DoFAccessor< structdim, DH >::DoFAccessor (  ) 

Default constructor. Provides an accessor that can't be used.

template<int structdim, class DH>
DoFAccessor< structdim, DH >::DoFAccessor ( const Triangulation< DH::dimension > *  tria,
const int  level,
const int  index,
const DH *  local_data 
)

Constructor


Member Function Documentation

template<int structdim, class DH>
void DoFAccessor< structdim, DH >::set_dof_handler ( DH *  dh  ) 

Reset the DoF handler pointer.

template<int structdim, class DH>
const DH& DoFAccessor< structdim, DH >::get_dof_handler (  )  const

Return a handle on the DoFHandler object which we are using.

template<int structdim, class DH>
DoFAccessor<structdim,DH>& DoFAccessor< structdim, DH >::operator= ( const DoFAccessor< structdim, DH > &  da  ) 

Copy operator.

template<int structdim, class DH>
void DoFAccessor< structdim, DH >::copy_from ( const DoFAccessor< structdim, DH > &  a  ) 

Implement the copy operator needed for the iterator classes.

template<int structdim, class DH>
TriaIterator<DH::dimension,DoFObjectAccessor<structdim,DH> > DoFAccessor< structdim, DH >::child ( const unsigned int  c  )  const

Return an iterator pointing to the the c-th child.

Reimplemented from TriaObjectAccessor< celldim, dim >.

template<int structdim, class DH>
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().

template<int structdim, class DH>
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

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().

template<int structdim, class DH>
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().

template<int structdim, class DH>
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

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().

template<int structdim, class DH>
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.

template<int structdim, class DH>
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.

template<int structdim, class DH>
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()).

template<int structdim, class DH>
const FiniteElement<DH::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.

template<int structdim, class DH>
bool DoFAccessor< structdim, DH >::operator== ( const DoFAccessor< structdim, DH > &   )  const [protected]

Compare for equality.

template<int structdim, class DH>
bool DoFAccessor< structdim, DH >::operator!= ( const DoFAccessor< structdim, DH > &   )  const [protected]

Compare for inequality.


Friends And Related Function Documentation

template<int structdim, class DH>
friend class TriaRawIterator [friend]

Iterator classes need to be friends because they need to access operator== and operator!=.

Reimplemented from TriaObjectAccessor< celldim, dim >.


Member Data Documentation

template<int structdim, class DH>
DH* DoFAccessor< structdim, DH >::dof_handler [protected]

Store the address of the DoFHandler object to be accessed.


The documentation for this class was generated from the following file:

deal.II documentation generated on Fri Nov 21 07:13:00 2008 by doxygen 1.5.6