Public Types | Public Member Functions | Static Public Attributes | Private Member Functions | Friends

DoFCellAccessor< DH > Class Template Reference
[Degrees of FreedomAccessor classes of the mesh iterators]

Inheritance diagram for DoFCellAccessor< DH >:
Inheritance graph
[legend]

List of all members.

Public Types

typedef DoFAccessor
< DH::dimension, DH >
::AccessorData 
AccessorData
typedef DoFAccessor
< DH::dimension, DH > 
BaseClass
typedef DH Container

Public Member Functions

internal::DoFHandler::Iterators
< DH >::cell_iterator 
parent () const
void set_dof_indices (const std::vector< unsigned int > &dof_indices)
void update_cell_dof_indices_cache () const
Constructors
 DoFCellAccessor (const Triangulation< DH::dimension, DH::space_dimension > *tria, const int level, const int index, const AccessorData *local_data)
template<int structdim2, int dim2, int spacedim2>
 DoFCellAccessor (const InvalidAccessor< structdim2, dim2, spacedim2 > &)
template<int dim2, class DH2 >
 DoFCellAccessor (const DoFAccessor< dim2, DH2 > &)
Accessing sub-objects and neighbors
internal::DoFHandler::Iterators
< DH >::cell_iterator 
neighbor (const unsigned int) const
internal::DoFHandler::Iterators
< DH >::cell_iterator 
child (const unsigned int) const
internal::DoFHandler::Iterators
< DH >::face_iterator 
face (const unsigned int i) const
internal::DoFHandler::Iterators
< DH >::cell_iterator 
neighbor_child_on_subface (const unsigned int face_no, const unsigned int subface_no) const
Extracting values from global vectors
template<class InputVector , typename number >
void get_dof_values (const InputVector &values, Vector< number > &local_values) const
template<class InputVector , typename ForwardIterator >
void get_dof_values (const InputVector &values, ForwardIterator local_values_begin, ForwardIterator local_values_end) const
template<class InputVector , typename ForwardIterator >
void get_dof_values (const ConstraintMatrix &constraints, const InputVector &values, ForwardIterator local_values_begin, ForwardIterator local_values_end) const
template<class OutputVector , typename number >
void set_dof_values (const Vector< number > &local_values, OutputVector &values) const
template<class InputVector , typename number >
void get_interpolated_dof_values (const InputVector &values, Vector< number > &interpolated_values) const
template<class OutputVector , typename number >
void set_dof_values_by_interpolation (const Vector< number > &local_values, OutputVector &values) const
template<typename number , typename OutputVector >
void distribute_local_to_global (const Vector< number > &local_source, OutputVector &global_destination) const
template<typename ForwardIterator , typename OutputVector >
void distribute_local_to_global (ForwardIterator local_source_begin, ForwardIterator local_source_end, OutputVector &global_destination) const
template<typename ForwardIterator , typename OutputVector >
void distribute_local_to_global (const ConstraintMatrix &constraints, ForwardIterator local_source_begin, ForwardIterator local_source_end, OutputVector &global_destination) const
template<typename number , typename OutputMatrix >
void distribute_local_to_global (const FullMatrix< number > &local_source, OutputMatrix &global_destination) const
template<typename number , typename OutputMatrix , typename OutputVector >
void distribute_local_to_global (const FullMatrix< number > &local_matrix, const Vector< number > &local_vector, OutputMatrix &global_matrix, OutputVector &global_vector) const
Accessing the DoF indices of this object
void get_dof_indices (std::vector< unsigned int > &dof_indices) const
Accessing the finite element associated with this object
const FiniteElement
< DH::dimension,
DH::space_dimension > & 
get_fe () const
unsigned int active_fe_index () const
void set_active_fe_index (const unsigned int i)

Static Public Attributes

static const unsigned int dim = DH::dimension
static const unsigned int spacedim = DH::space_dimension

Private Member Functions

DoFCellAccessor< DH > & operator= (const DoFCellAccessor< DH > &da)

Friends

class DoFHandler
struct internal::DoFCellAccessor::Implementation

Detailed Description

template<class DH>
class DoFCellAccessor< DH >

Grant access to the degrees of freedom on a cell.

Note that since for the class we derive from, i.e. DoFAccessor<dim>, the two template parameters are equal, the base class is actually derived from CellAccessor, which makes the functions of this class available to the DoFCellAccessor class as well.

Author:
Wolfgang Bangerth, 1998

Definition at line 1398 of file dof_accessor.h.


Member Typedef Documentation

template<class DH>
typedef DoFAccessor<DH::dimension,DH>::AccessorData DoFCellAccessor< DH >::AccessorData

Declare the data type that this accessor class expects to get passed from the iterator classes.

Reimplemented from DoFAccessor< DH::dimension, DH >.

Definition at line 1417 of file dof_accessor.h.

template<class DH>
typedef DoFAccessor<DH::dimension,DH> DoFCellAccessor< DH >::BaseClass

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

Reimplemented from DoFAccessor< DH::dimension, DH >.

Definition at line 1425 of file dof_accessor.h.

template<class DH>
typedef DH DoFCellAccessor< DH >::Container

Define the type of the container this is part of.

Definition at line 1431 of file dof_accessor.h.


Constructor & Destructor Documentation

template<class DH>
DoFCellAccessor< DH >::DoFCellAccessor ( const Triangulation< DH::dimension, DH::space_dimension > *  tria,
const int  level,
const int  index,
const AccessorData local_data 
)

Constructor

template<class DH>
template<int structdim2, int dim2, int spacedim2>
DoFCellAccessor< DH >::DoFCellAccessor ( 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).

template<class DH>
template<int dim2, class DH2 >
DoFCellAccessor< DH >::DoFCellAccessor ( const DoFAccessor< dim2, DH2 > &   )

Another conversion operator between objects that don't make sense, just like the previous one.


Member Function Documentation

template<class DH>
internal::DoFHandler::Iterators<DH>::cell_iterator DoFCellAccessor< DH >::parent (  ) const

Return the parent as a DoF cell iterator. This function is needed since the parent function of the base class returns a cell accessor without access to the DoF data.

Reimplemented from DoFAccessor< DH::dimension, DH >.

template<class DH>
internal::DoFHandler::Iterators<DH>::cell_iterator DoFCellAccessor< DH >::neighbor ( const unsigned  int ) const

Return the ith neighbor as a DoF cell iterator. This function is needed since the neighbor function of the base class returns a cell accessor without access to the DoF data.

template<class DH>
internal::DoFHandler::Iterators<DH>::cell_iterator DoFCellAccessor< DH >::child ( const unsigned  int ) const

Return the ith child as a DoF cell iterator. This function is needed since the child function of the base class returns a cell accessor without access to the DoF data.

Reimplemented from DoFAccessor< DH::dimension, DH >.

template<class DH>
internal::DoFHandler::Iterators<DH>::face_iterator DoFCellAccessor< DH >::face ( const unsigned int  i ) const

Return an iterator to the ith face of this cell.

This function is not implemented in 1D, and maps to DoFAccessor::line in 2D.

template<class DH>
internal::DoFHandler::Iterators<DH>::cell_iterator DoFCellAccessor< DH >::neighbor_child_on_subface ( const unsigned int  face_no,
const unsigned int  subface_no 
) const

Return the result of the neighbor_child_on_subface function of the base class, but convert it so that one can also access the DoF data (the function in the base class only returns an iterator with access to the triangulation data).

template<class DH>
template<class InputVector , typename number >
void DoFCellAccessor< DH >::get_dof_values ( const InputVector &  values,
Vector< number > &  local_values 
) const

Return the values of the given vector restricted to the dofs of this cell 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 only callable for active cells.

The input vector may be either a Vector<float>, Vector<double>, or a BlockVector<double>, or a PETSc or Trilinos vector if deal.II is compiled to support these libraries. It is in the responsibility of the caller to assure that the types of the numbers stored in input and output vectors are compatible and with similar accuracy.

template<class DH>
template<class InputVector , typename ForwardIterator >
void DoFCellAccessor< DH >::get_dof_values ( const InputVector &  values,
ForwardIterator  local_values_begin,
ForwardIterator  local_values_end 
) const

Return the values of the given vector restricted to the dofs of this cell 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 only callable for active cells.

The input vector may be either a Vector<float>, Vector<double>, or a BlockVector<double>, or a PETSc or Trilinos vector if deal.II is compiled to support these libraries. It is in the responsibility of the caller to assure that the types of the numbers stored in input and output vectors are compatible and with similar accuracy.

template<class DH>
template<class InputVector , typename ForwardIterator >
void DoFCellAccessor< DH >::get_dof_values ( const ConstraintMatrix constraints,
const InputVector &  values,
ForwardIterator  local_values_begin,
ForwardIterator  local_values_end 
) const

Return the values of the given vector restricted to the dofs of this cell 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 only callable for active cells.

The input vector may be either a Vector<float>, Vector<double>, or a BlockVector<double>, or a PETSc or Trilinos vector if deal.II is compiled to support these libraries. It is in the responsibility of the caller to assure that the types of the numbers stored in input and output vectors are compatible and with similar accuracy. The ConstraintMatrix passed as an argument to this function makes sure that constraints are correctly distributed when the dof values are calculated.

template<class DH>
template<class OutputVector , typename number >
void DoFCellAccessor< DH >::set_dof_values ( const Vector< number > &  local_values,
OutputVector &  values 
) const

This function is the counterpart to get_dof_values(): it takes a vector of values for the degrees of freedom of the cell pointed to by this iterator and writes these values into the global data vector values. This function is only callable for active cells.

Note that for continuous finite elements, calling this function affects the dof values on neighboring cells as well. It may also violate continuity requirements for hanging nodes, if neighboring cells are less refined than the present one. These requirements are not taken care of and must be enforced by the user afterwards.

The vector has to have the right size before being passed to this function.

The output vector may be either a Vector<float>, Vector<double>, or a BlockVector<double>, or a PETSc vector if deal.II is compiled to support these libraries. It is in the responsibility of the caller to assure that the types of the numbers stored in input and output vectors are compatible and with similar accuracy.

template<class DH>
template<class InputVector , typename number >
void DoFCellAccessor< DH >::get_interpolated_dof_values ( const InputVector &  values,
Vector< number > &  interpolated_values 
) const

Return the interpolation of the given finite element function to the present cell. In the simplest case, the cell is a terminal one, i.e. has no children; then, the returned value is the vector of nodal values on that cell. You could then as well get the desired values through the get_dof_values function. In the other case, when the cell has children, we use the restriction matrices provided by the finite element class to compute the interpolation from the children to the present cell.

It is assumed that both vectors already have the right size beforehand.

Unlike the get_dof_values() function, this function works on cells rather than to lines, quads, and hexes, since interpolation is presently only provided for cells by the finite element classes.

template<class DH>
template<class OutputVector , typename number >
void DoFCellAccessor< DH >::set_dof_values_by_interpolation ( const Vector< number > &  local_values,
OutputVector &  values 
) const

This, again, is the counterpart to get_interpolated_dof_values(): you specify the dof values on a cell and these are interpolated to the children of the present cell and set on the terminal cells.

In principle, it works as follows: if the cell pointed to by this object is terminal, then the dof values are set in the global data vector by calling the set_dof_values() function; otherwise, the values are prolonged to each of the children and this function is called for each of them.

Using the get_interpolated_dof_values() and this function, you can compute the interpolation of a finite element function to a coarser grid by first getting the interpolated solution on a cell of the coarse grid and afterwards redistributing it using this function.

Note that for continuous finite elements, calling this function affects the dof values on neighboring cells as well. It may also violate continuity requirements for hanging nodes, if neighboring cells are less refined than the present one, or if their children are less refined than the children of this cell. These requirements are not taken care of and must be enforced by the user afterward.

It is assumed that both vectors already have the right size beforehand. This function relies on the existence of a natural interpolation property of finite element spaces of a cell to its children, denoted by the prolongation matrices of finite element classes. For some elements, the spaces on coarse and fine grids are not nested, in which case the interpolation to a child is not the identity; refer to the documentation of the respective finite element class for a description of what the prolongation matrices represent in this case.

Unlike the set_dof_values() function, this function is associated to cells rather than to lines, quads, and hexes, since interpolation is presently only provided for cells by the finite element objects.

The output vector may be either a Vector<float>, Vector<double>, or a BlockVector<double>, or a PETSc vector if deal.II is compiled to support these libraries. It is in the responsibility of the caller to assure that the types of the numbers stored in input and output vectors are compatible and with similar accuracy.

template<class DH>
template<typename number , typename OutputVector >
void DoFCellAccessor< DH >::distribute_local_to_global ( const Vector< number > &  local_source,
OutputVector &  global_destination 
) const

Distribute a local (cell based) vector to a global one by mapping the local numbering of the degrees of freedom to the global one and entering the local values into the global vector.

The elements are added up to the elements in the global vector, rather than just set, since this is usually what one wants.

template<class DH>
template<typename ForwardIterator , typename OutputVector >
void DoFCellAccessor< DH >::distribute_local_to_global ( ForwardIterator  local_source_begin,
ForwardIterator  local_source_end,
OutputVector &  global_destination 
) const

Distribute a local (cell based) vector in iterator format to a global one by mapping the local numbering of the degrees of freedom to the global one and entering the local values into the global vector.

The elements are added up to the elements in the global vector, rather than just set, since this is usually what one wants.

template<class DH>
template<typename ForwardIterator , typename OutputVector >
void DoFCellAccessor< DH >::distribute_local_to_global ( const ConstraintMatrix constraints,
ForwardIterator  local_source_begin,
ForwardIterator  local_source_end,
OutputVector &  global_destination 
) const

Distribute a local (cell based) vector in iterator format to a global one by mapping the local numbering of the degrees of freedom to the global one and entering the local values into the global vector.

The elements are added up to the elements in the global vector, rather than just set, since this is usually what one wants. Moreover, the ConstraintMatrix passed to this function makes sure that also constraints are eliminated in this process.

template<class DH>
template<typename number , typename OutputMatrix >
void DoFCellAccessor< DH >::distribute_local_to_global ( const FullMatrix< number > &  local_source,
OutputMatrix &  global_destination 
) const

This function does much the same as the distribute_local_to_global(Vector,Vector) function, but operates on matrices instead of vectors. If the matrix type is a sparse matrix then it is supposed to have non-zero entry slots where required.

template<class DH>
template<typename number , typename OutputMatrix , typename OutputVector >
void DoFCellAccessor< DH >::distribute_local_to_global ( const FullMatrix< number > &  local_matrix,
const Vector< number > &  local_vector,
OutputMatrix &  global_matrix,
OutputVector &  global_vector 
) const

This function does what the two distribute_local_to_global functions with vector and matrix argument do, but all at once.

template<class DH>
void DoFCellAccessor< DH >::get_dof_indices ( std::vector< unsigned int > &  dof_indices ) const

Return the indices of the dofs of this quad 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.

It is assumed that the vector already has the right size beforehand.

This function reimplements the same function in the base class. The functions in the base classes are available for all geometric objects, i.e. even in 3d they can be used to access the dof indices of edges, for example. On the other hand, the most common case is clearly the use on cells, which is why we cache the array for each cell, but not edge. To retrieve the cached values, rather than collect the necessary information every time, this function overwrites the one in the base class.

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.

template<class DH>
const FiniteElement<DH::dimension,DH::space_dimension>& DoFCellAccessor< DH >::get_fe (  ) const

Return the finite element that is used on the cell pointed to by this iterator. For non-hp DoF handlers, this is of course always the same element, independent of the cell we are presently on, but for hp DoF handlers, this may change from cell to cell.

template<class DH>
unsigned int DoFCellAccessor< DH >::active_fe_index (  ) const

Returns the index inside the hp::FECollection of the FiniteElement used for this cell.

template<class DH>
void DoFCellAccessor< DH >::set_active_fe_index ( const unsigned int  i )

Sets the index of the FiniteElement used for this cell.

template<class DH>
void DoFCellAccessor< DH >::set_dof_indices ( const std::vector< unsigned int > &  dof_indices )

Set the DoF indices of this cell to the given values. This function bypasses the DoF cache, if one exists for the given DoF handler class.

template<class DH>
void DoFCellAccessor< DH >::update_cell_dof_indices_cache (  ) const

Update the cache in which we store the dof indices of this cell.

template<class DH>
DoFCellAccessor<DH>& DoFCellAccessor< DH >::operator= ( const DoFCellAccessor< DH > &  da ) [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.


Friends And Related Function Documentation

template<class DH>
friend class DoFHandler [friend]

Make the DoFHandler class a friend so that it can call the update_cell_dof_indices_cache() function

Reimplemented from DoFAccessor< DH::dimension, DH >.

Definition at line 2070 of file dof_accessor.h.

template<class DH>
friend struct internal::DoFCellAccessor::Implementation [friend]

Reimplemented from DoFAccessor< DH::dimension, DH >.

Definition at line 2071 of file dof_accessor.h.


Member Data Documentation

template<class DH>
const unsigned int DoFCellAccessor< DH >::dim = DH::dimension [static]

Extract dimension from DH.

Definition at line 1404 of file dof_accessor.h.

template<class DH>
const unsigned int DoFCellAccessor< DH >::spacedim = DH::space_dimension [static]

Extract space dimension from DH.

Definition at line 1409 of file dof_accessor.h.


The documentation for this class was generated from the following file:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

deal.II documentation generated on Fri Feb 3 2012 06:04:07 by doxygen 1.7.2