#define __deal2__dof_accessor_h
-/*------------------------- dof_iterator.h ------------------------*/
-
-
#include <base/config.h>
#include <grid/tria_accessor.h>
#include <vector>
/**
* Define the basis for accessors to the degrees of freedom.
*
- * Note that it is allowed to construct an object of which the
- * @p dof_handler pointer is a Null pointer. Such an object would
+ * Note that it is allowed to construct an object of which
+ * #dof_handler is a Null pointer. Such an object would
* result in a strange kind of behaviour, though every reasonable
* operating system should disallow access through that pointer.
* The reason we do not check for the null pointer in the
* which has an invalid dof handler pointer. This is to guarantee
* that every iterator which is once assigned to is a valid
* object. However, this assertion only holds in debug mode, when
- * the @p Assert macro is switched on.
+ * the #Assert macro is switched on.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998
*/
template <int dim>
/**
* This should be the default constructor.
- * We cast away the @p constness of the
+ * We cast away the constness of the
* pointer which clearly is EVIL but
* we can't help without making all
* functions which could somehow use
/**
* Exception for child classes
+ *
+ * @ingroup Exceptions
*/
DeclException0 (ExcInvalidObject);
/**
* Exception
+ *
+ * @ingroup Exceptions
*/
DeclException0 (ExcVectorNotEmpty);
/**
* Exception
+ *
+ * @ingroup Exceptions
*/
DeclException0 (ExcVectorDoesNotMatch);
/**
* Exception
+ *
+ * @ingroup Exceptions
*/
DeclException0 (ExcMatrixDoesNotMatch);
/**
- * Exception
+ * A function has been called for
+ * a cell which should be active,
+ * but is refined. @ref GlossActive
+ *
+ * @ingroup Exceptions
*/
DeclException0 (ExcNotActive);
* the inheritance is automatically chosen to be from CellAccessor if the
* object under consideration has full dimension, i.e. constitutes a cell.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1999
*/
template <int celldim, int dim>
* the inheritance is automatically chosen to be from CellAccessor if the
* object under consideration has full dimension, i.e. constitutes a cell.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1999
*/
template <int dim>
* lines in 1D-, 2D-, etc dimensions).
*
*
- * @sect3{Usage}
+ * <h3>Usage</h3>
*
* The DoFDimensionInfo classes inherited by the DoFHandler classes
* declare typedefs to iterators using the accessors declared in this class
* as they provide easier typing (much less complicated names!).
*
*
- * @sect3{Notes about the class hierarchy structure}
+ * <h3>Notes about the class hierarchy structure</h3>
*
* See the report on this subject, which is available from the general
* documentation directory.
*
- *
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998; Guido Kanschat, 1999
*
* (Internal: inheritance is necessary for the general template due to
const AccessorData *local_data);
/**
- * Index of the @p ith degree
+ * Index of the <i>i</i>th degree
* of freedom of this object.
*/
unsigned int dof_index (const unsigned int i) const;
/**
- * Set the index of the @p ith degree
+ * Set the index of the <i>i</i>th degree
* of freedom of this object to @p index.
*/
void set_dof_index (const unsigned int i,
const int index) const;
/**
- * Index of the @p ith degree
+ * Index of the <i>i</i> degree
* on the @p vertexth vertex.
*/
unsigned int vertex_dof_index (const unsigned int vertex,
const unsigned int i) const;
/**
- * Set the index of the @p ith degree
+ * Set the index of the <i>i</i> degree
* on the @p vertexth vertex to @p index.
*/
void set_vertex_dof_index (const unsigned int vertex,
/**
* This function is the counterpart to
- * @p get_dof_values: it takes a vector
+ * 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
OutputVector &values) const;
/**
- * Pointer to the @p ith line
+ * Pointer to the <i>i</i>th line
* bounding this Object.
*/
TriaIterator<dim,DoFObjectAccessor<1, dim> >
line (const unsigned int i) const;
/**
- * Pointer to the @p ith quad
+ * Pointer to the <i>i</i>th quad
* bounding this Object.
*/
TriaIterator<dim,DoFObjectAccessor<2, dim> >
quad (const unsigned int i) const;
/**
- * @p ith child as a DoFObjectAccessor
+ * <i>i</i>th child as a DoFObjectAccessor
* iterator. This function is needed since
* the child function of the base
* class returns a hex accessor without
* lines in 1D-, 2D-, etc dimensions).
*
*
- * @sect3{Usage}
+ * <h3>Usage</h3>
*
* The DoFDimensionInfo classes inherited by the DoFHandler classes
* declare typedefs to iterators using the accessors declared in this class
* as they provide easier typing (much less complicated names!).
*
*
- * @sect3{Notes about the class hierarchy structure}
+ * <h3>Notes about the class hierarchy structure</h3>
*
* Inheritance from <tt>DoFObjectAccessor_Inheritance<1,dim>::BaseClass</tt> yields
* inheritance from <tt>CellAccessor<1></tt> if <tt>dim==1</tt> and from
* of this class shares all features of cells in one dimension, but behaves
* like an ordinary line in all other cases.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998
*/
template <int dim>
/**
* This function is the counterpart to
- * @p get_dof_values: it takes a vector
+ * 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
/**
* Grant access to the degrees of freedom located on quads.
*
- * @ref DoFObjectAccessor
+ * @ingroup dofs
*/
template <int dim>
class DoFObjectAccessor<2, dim> : public DoFAccessor<dim>,
/**
* This function is the counterpart to
- * @p get_dof_values: it takes a vector
+ * 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
/**
* Grant access to the degrees of freedom located on hexes.
*
- * @ref DoFObjectAccessor
+ * @ingroup dofs
*/
template <int dim>
class DoFObjectAccessor<3, dim> : public DoFAccessor<dim>,
/**
* This function is the counterpart to
- * @p get_dof_values: it takes a vector
+ * 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
* CellAccessor, which makes the functions of this class available to the
* DoFCellAccessor class as well.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998
*/
template <int dim>
* decided what the this function
* does in these cases.
*
- * Unlike the @p get_dof_values
+ * Unlike the get_dof_values()
* function, this function is
* associated to cells rather
* than to lines, quads, and
/**
* This, again, is the
* counterpart to
- * @p get_interpolated_dof_values:
+ * get_interpolated_dof_values():
* you specify the dof values on
* a cell and these are
* interpolated to the children
* to by this object is terminal,
* then the dof values are set in
* the global data vector by
- * calling the @p set_dof_values
+ * calling the set_dof_values()
* function; otherwise, the
* values are prolonged to each
* of the children and this
* them.
*
* Using the
- * @p get_interpolated_dof_values
+ * get_interpolated_dof_values()
* and this function, you can
* compute the interpolation of a
* finite element function to a
* what the prolongation matrices
* represent in this case.
*
- * Unlike the @p set_dof_values
+ * Unlike the set_dof_values()
* function, this function is
* associated to cells rather
* than to lines, quads, and
* which has to applied to the solution vector.
*
*
- * @sect3{Condensing matrices and sparsity patterns}
+ * <h3>Condensing matrices and sparsity patterns</h3>
*
* Condensation of a matrix is done in four steps: first one builds the
* sparsity pattern (e.g. using
- * DoFHandler@p ::create_sparsity_pattern); then the sparsity pattern
+ * DoFHandler::create_sparsity_pattern); then the sparsity pattern
* of the condensed matrix is made out of the original sparsity pattern and
* the constraints; third, the global matrix is assembled; and fourth, the
* matrix is finally condensed. To do these steps, you have (at least) two
* in the already condensed form.
*
*
- * @sect3{Condensing vectors}
+ * <h3>Condensing vectors</h3>
*
* Condensing vectors works exactly as described above for matrices. Note that
* condensation is an idempotent operation, i.e. doing it more than once on a
* techniques as mentioned above to avoid their use.
*
*
- * @sect3{Avoiding explicit condensation}
+ * <h3>Avoiding explicit condensation</h3>
*
* Sometimes, one wants to avoid condensation at all. This may be the case
* since condensation is an expensive operation, or because no condense()
* system.
*
*
- * @sect3{Distributing constraints}
+ * <h3>Distributing constraints</h3>
*
* After solving the condensed system of equations, the solution vector has to
* be redistributed. This is done by the two @p distribute function, one
* for unconstrained nodes, and constrained nodes need to get their values in
* a second step.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998, 2004
*/
class ConstraintMatrix : public Subscriptor
/**
* Exception
+ *
+ * @ingroup Exceptions
*/
DeclException0 (ExcMatrixIsClosed);
/**
* Exception
+ *
+ * @ingroup Exceptions
*/
DeclException0 (ExcMatrixNotClosed);
/**
* Exception
+ *
+ * @ingroup Exceptions
*/
DeclException1 (ExcLineInexistant,
unsigned int,
<< " does not exist.");
/**
* Exception
+ *
+ * @ingroup Exceptions
*/
DeclException0 (ExcWrongDimension);
/**
* Exception
+ *
+ * @ingroup Exceptions
*/
DeclException4 (ExcEntryAlreadyExists,
int, int, double, double,
<< "by " << (arg4-arg3) << ".");
/**
* Exception
+ *
+ * @ingroup Exceptions
*/
DeclException2 (ExcDoFConstrainedToConstrainedDoF,
int, int,
<< ", but that one is also constrained. This is not allowed!");
/**
* Exception.
+ *
+ * @ingroup Exceptions
*/
DeclException1 (ExcDoFIsConstrainedFromBothObjects,
int,
<< " is constrained from both object in a merge operation.");
/**
* Exception
+ *
+ * @ingroup Exceptions
*/
DeclException1 (ExcDoFIsConstrainedToConstrainedDoF,
int,
* @ref DoFDimensionInfo<1>
* @ref DoFDimensionInfo<2>
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998
*/
template <int dim>
class DoFDimensionInfo
{
+ public:
+ /**
+ * For internal use only
+ */
+ typedef void * raw_line_iterator;
+ /**
+ * Iterator for lines in the
+ * DoFHandler
+ */
+ typedef void * line_iterator;
+ /**
+ * Iterator skipping all refined
+ * lines
+ */
+ typedef void * active_line_iterator;
+
+ /**
+ * For internal use only
+ */
+ typedef void * raw_quad_iterator;
+ /**
+ * Iterator for quadrilaterals in
+ * the DoFHandler
+ */
+ typedef void * quad_iterator;
+ /**
+ * Iterator skipping all refined
+ * quadrilaterals
+ */
+ typedef void * active_quad_iterator;
+
+ /**
+ * For internal use only
+ */
+ typedef void * raw_hex_iterator;
+ /**
+ * Iterator for hexahedra in the
+ * DoFHandler
+ */
+ typedef void * hex_iterator;
+ /**
+ * Iterator skipping all refined
+ * hexahedra
+ */
+ typedef void * active_hex_iterator;
+
+ /**
+ * For internal use only
+ */
+ typedef void * raw_cell_iterator;
+ /**
+ * All cells of the DoFHandler;
+ * depending on the dimension
+ * this is typedefed to
+ * <ol>
+ * <li> #line_iterator
+ * <li> #quad_iterator
+ * <li> #hex_iterator
+ * </ol>
+ */
+ typedef void * cell_iterator;
+ /**
+ * All active cells of the DoFHandler
+ * Depending on the dimension
+ * this is typedefed to
+ * <ol>
+ * <li> #active_line_iterator
+ * <li> #active_quad_iterator
+ * <li> #active_hex_iterator
+ * </ol>
+ */
+ typedef void * active_cell_iterator;
+
+ /**
+ * For internal use only
+ */
+ typedef void * raw_face_iterator;
+ /**
+ * All faces of the DoFHandler
+ * Depending on the dimension
+ * this is typedefed to
+ * <ol>
+ * <li> <tt>void*</tt>
+ * <li> #line_iterator
+ * <li> #quad_iterator
+ * </ol>
+ */
+ typedef void * face_iterator;
+ /**
+ * All active faces of the DoFHandler
+ * Depending on the dimension
+ * this is typedefed to
+ * <ol>
+ * <li> <tt>void*</tt>
+ * <li> #active_line_iterator
+ * <li> #active_quad_iterator
+ * </ol>
+ */
+ typedef void * active_face_iterator;
};
*
* The types have the same meaning as those declared in TriaDimensionInfo<2>.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998
*/
template <>
*
* The types have the same meaning as those declared in TriaDimensionInfo<2>.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998
*/
template <>
*
* The types have the same meaning as those declared in TriaDimensionInfo<3>.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998
*/
template <>
* global matrices as well as matrices living on (parts of) the boundary.
*
*
- * @sect3{Distribution of indices for degrees of freedom}
+ * <h3>Distribution of indices for degrees of freedom</h3>
*
* The degrees of freedom (`dofs') are distributed on the given triangulation
- * by the function <tt>distribute_dofs()</tt>. It gets passed a finite element object
+ * by the function distribute_dofs(). It gets passed a finite element object
* describing how many degrees of freedom are located on vertices, lines, etc.
* It traverses the triangulation cell by cell and numbers the dofs of that
* cell if not yet numbered. For non-multigrid algorithms, only active cells
* the implemented algorithms.
*
*
- * @sect3{User defined renumbering schemes}
+ * <h3>User defined renumbering schemes</h3>
*
* The DoFRenumbering class offers a number of renumbering
* schemes like the Cuthill-McKey scheme. Basically, the function sets
* up an array in which for each degree of freedom the index is stored
* which is to be assigned by the renumbering. Using this array, the
- * <tt>renumber_dofs(vector<unsigned int>)</tt> function is called, which actually
+ * renumber_dofs() function is called, which actually
* does the change from old DoF indices to the ones given in the
* array. In some cases, however, a user may want to compute her own
* renumbering order; in this case, allocate an array with one element
* <tt>renumber_dofs(vector<unsigned int>)</tt> with the array, which converts old
* into new degree of freedom indices.
*
- *
- *
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998
*/
template <int dim>
*/
std::vector<unsigned int> vertex_dofs;
- /**
+ /*
* Make accessor objects friends.
*/
template <int dim1> friend class DoFAccessor;
- /**
+ /*
* Make accessor objects friends.
*/
template <int dim1, int dim2> friend class DoFObjectAccessor;
* class, but do not actually use it. Rather, only specializations of
* this class are used.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998
*/
template <int N>
* Store the indices of the degrees of freedom which are located on
* the lines.
*
- * @sect3{Information for all DoFLevel classes}
+ * <h3>Information for all DoFLevel classes</h3>
*
* The <tt>DoFLevel<N></tt> classes
* store the global indices of the degrees of freedom for each cell on a
* is used, the indices are stored in the @p vertex_dofs array of the
* DoFHandler class.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998
*/
template <>
* Store the indices of the degrees of freedom which are located on
* quads. See DoFLevel<1> for more information.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998
*/
template <>
* Store the indices of the degrees of freedom which are located on
* hexhedra. See DoFLevel<1> for more information.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, 1998
*/
template <>
* Implementation of a number of renumbering algorithms for the degrees of
* freedom on a triangulation.
*
- * @sect2{Cuthill-McKee like algorithms}
- *
- * Within this class, the Cuthill-McKee algorithm is implemented. It starts
- * at a degree of freedom, searches the other DoFs for those which are couple
- * with the one we started with and numbers these in a certain way. It then
- * finds the second level of DoFs, namely those that couple with those of
- * the previous level (which were those that coupled with the initial DoF)
- * and numbers these. And so on. For the details of the algorithm, especially
- * the numbering within each level, we refer the reader to the book of
- * Schwarz (H.R.Schwarz: Methode der finiten Elemente). The reverse Cuthill-McKee
- * algorithm does the same job, but numbers all elements in the reverse order.
- *
- * These algorithms
- * have one major drawback: they require a good starting point, i.e. the degree
- * of freedom index afterwards to be numbered zero. This can thus be given by
- * the user, e.g. by exploiting knowledge of the actual topology of the
- * domain. It is also possible to give several starting indices, which may
- * be used to simulate a simple upstream numbering (by giving the inflow
- * dofs as starting values) or to make preconditioning faster (by letting
+ * <h3>Cuthill-McKee like algorithms</h3>
+ *
+ * Within this class, the Cuthill-McKee algorithm is implemented. It
+ * starts at a degree of freedom, searches the other DoFs for those
+ * which are couple with the one we started with and numbers these in
+ * a certain way. It then finds the second level of DoFs, namely those
+ * that couple with those of the previous level (which were those that
+ * coupled with the initial DoF) and numbers these. And so on. For the
+ * details of the algorithm, especially the numbering within each
+ * level, we refer the reader to the book of Schwarz (H.R.Schwarz:
+ * Methode der finiten Elemente). The reverse Cuthill-McKee algorithm
+ * does the same job, but numbers all elements in the reverse order.
+ *
+ * These algorithms have one major drawback: they require a good
+ * starting point, i.e. the degree of freedom index afterwards to be
+ * numbered zero. This can thus be given by the user, e.g. by
+ * exploiting knowledge of the actual topology of the domain. It is
+ * also possible to give several starting indices, which may be used
+ * to simulate a simple upstream numbering (by giving the inflow dofs
+ * as starting values) or to make preconditioning faster (by letting
* the dirichlet boundary indices be starting points).
*
- * If no starting index is given, one is chosen by the program, namely one
- * with the smallest coordination number (the coordination number is the
- * number of other dofs this dof couples with). This dof is usually located
- * on the boundary of the domain. There is, however, large ambiguity in this
- * when using the hierarchical meshes used in this library, since in most
- * cases the computational domain is not approximated by tilting and deforming
- * elements and by plugging together variable numbers of elements at vertices,
- * but rather by hierarchical refinement. There is therefore a large number
- * of dofs with equal coordination numbers. The renumbering algorithms will
+ * If no starting index is given, one is chosen by the program, namely
+ * one with the smallest coordination number (the coordination number
+ * is the number of other dofs this dof couples with). This dof is
+ * usually located on the boundary of the domain. There is, however,
+ * large ambiguity in this when using the hierarchical meshes used in
+ * this library, since in most cases the computational domain is not
+ * approximated by tilting and deforming elements and by plugging
+ * together variable numbers of elements at vertices, but rather by
+ * hierarchical refinement. There is therefore a large number of dofs
+ * with equal coordination numbers. The renumbering algorithms will
* therefore not give optimal results.
*
- * In the book of Schwarz (H.R.Schwarz: Methode der finiten Elemente), it is
- * advised to test many starting points, if possible all with the smallest
- * coordination number and also those with slightly higher numbers. However,
- * this seems only possible for meshes with at most several dozen or a few
- * hundred elements found in small engineering problems of the early 1980s
- * (the second edition was published in 1984), but certainly not with those
- * used in this library, featuring several 10,000 to a few 100,000 elements.
- *
- * On the other hand, the need to reduce the bandwidth has decreased since
- * with the mentioned number of cells, only iterative solution methods are
- * able to solve the resulting matrix systems. These, however, are not so
- * demanding with respect to the bandwidth as direct solvers used for
- * smaller problems. Things like upstream numbering become much more important
- * in recent times, so the suboptimality of the renumbering algorithms is
- * not that important any more.
+ * In the book of Schwarz (H.R.Schwarz: Methode der finiten Elemente),
+ * it is advised to test many starting points, if possible all with
+ * the smallest coordination number and also those with slightly
+ * higher numbers. However, this seems only possible for meshes with
+ * at most several dozen or a few hundred elements found in small
+ * engineering problems of the early 1980s (the second edition was
+ * published in 1984), but certainly not with those used in this
+ * library, featuring several 10,000 to a few 100,000 elements.
+ *
+ * On the other hand, the need to reduce the bandwidth has decreased
+ * since with the mentioned number of cells, only iterative solution
+ * methods are able to solve the resulting matrix systems. These,
+ * however, are not so demanding with respect to the bandwidth as
+ * direct solvers used for smaller problems. Things like upstream
+ * numbering become much more important in recent times, so the
+ * suboptimality of the renumbering algorithms is not that important
+ * any more.
*
*
- * @sect3{Implementation of renumbering schemes}
- *
- * The renumbering algorithms need quite a lot of memory, since they have
- * to store for each dof with which other dofs it couples. This is done
- * using a SparsityPattern object used to store the sparsity pattern of
- * matrices. It
- * is not useful for the user to do anything between distributing the dofs
- * and renumbering, i.e. the calls to DoFHandler@p ::distribute_dofs and
- * DoFHandler@p ::renumber_dofs should follow each other immediately. If
- * you try to create a sparsity pattern or anything else in between, these
- * will be invalid afterwards.
- *
- * The renumbering may take care of dof-to-dof couplings only induced by
- * eliminating constraints. In addition to the memory consumption mentioned
- * above, this also takes quite some computational time, but it may be
- * switched off upon calling the @p renumber_dofs function. This will then
- * give inferior results, since knots in the graph (representing dofs)
- * are not found to be neighbors even if they would be after condensation.
+ * <h4>Implementation of renumbering schemes</h4>
+ *
+ * The renumbering algorithms need quite a lot of memory, since they
+ * have to store for each dof with which other dofs it couples. This
+ * is done using a SparsityPattern object used to store the sparsity
+ * pattern of matrices. It is not useful for the user to do anything
+ * between distributing the dofs and renumbering, i.e. the calls to
+ * DoFHandler::distribute_dofs and DoFHandler::renumber_dofs should
+ * follow each other immediately. If you try to create a sparsity
+ * pattern or anything else in between, these will be invalid
+ * afterwards.
+ *
+ * The renumbering may take care of dof-to-dof couplings only induced
+ * by eliminating constraints. In addition to the memory consumption
+ * mentioned above, this also takes quite some computational time, but
+ * it may be switched off upon calling the @p renumber_dofs
+ * function. This will then give inferior results, since knots in the
+ * graph (representing dofs) are not found to be neighbors even if
+ * they would be after condensation.
*
- * The renumbering algorithms work on a purely algebraic basis, due to the
- * isomorphism between the graph theoretical groundwork underlying the
- * algorithms and binary matrices (matrices of which the entries are binary
- * values) represented by the sparsity patterns. In special, the algorithms
- * do not try to exploit topological knowledge (e.g. corner detection) to
- * find appropriate starting points. This way, however, they work in
- * arbitrary space dimension.
- *
- * If you want to give starting points, you may give a list of dof indices
- * which will form the first step of the renumbering. The dofs of the list
- * will be consecutively numbered starting with zero, i.e. this list is not
- * renumbered according to the coordination number of the nodes. Indices not
- * in the allowed range are deleted. If no index is allowed, the algorithm
- * will search for its own starting point.
+ * The renumbering algorithms work on a purely algebraic basis, due to
+ * the isomorphism between the graph theoretical groundwork underlying
+ * the algorithms and binary matrices (matrices of which the entries
+ * are binary values) represented by the sparsity patterns. In
+ * special, the algorithms do not try to exploit topological knowledge
+ * (e.g. corner detection) to find appropriate starting points. This
+ * way, however, they work in arbitrary space dimension.
+ *
+ * If you want to give starting points, you may give a list of dof
+ * indices which will form the first step of the renumbering. The dofs
+ * of the list will be consecutively numbered starting with zero,
+ * i.e. this list is not renumbered according to the coordination
+ * number of the nodes. Indices not in the allowed range are
+ * deleted. If no index is allowed, the algorithm will search for its
+ * own starting point.
*
*
- * @sect3{Results of renumbering}
- *
- * The renumbering schemes mentioned above do not lead to optimal results.
- * However, after all there is no algorithm that accomplishes this within
- * reasonable time. There are situations where the lack of optimality even
- * leads to worse results than with the original, crude, levelwise numering
- * scheme; one of these examples is a mesh of four cells of which always
- * those cells are refined which are neighbors to the center (you may call
- * this mesh a `zoom in' mesh). In one such example the bandwidth was
- * increased by about 50 per cent.
+ * <h4>Results of renumbering</h4>
+ *
+ * The renumbering schemes mentioned above do not lead to optimal
+ * results. However, after all there is no algorithm that
+ * accomplishes this within reasonable time. There are situations
+ * where the lack of optimality even leads to worse results than with
+ * the original, crude, levelwise numering scheme; one of these
+ * examples is a mesh of four cells of which always those cells are
+ * refined which are neighbors to the center (you may call this mesh a
+ * `zoom in' mesh). In one such example the bandwidth was increased by
+ * about 50 per cent.
*
* In most other cases, the bandwith is reduced significantly. The reduction
* is the better the less structured the grid is. With one grid where the
* may be difficult, however, and in many cases will not justify the effort.
*
*
- * @sect2{Component-wise numbering}
+ * <h3>Component-wise numbering</h3>
*
* For finite elements composed of several base elements using the FESystem
* class, or for elements which provide several components themselves, it
* different components.
*
* This kind of numbering may be obtained by calling the
- * @p component_wise function of this class. Since it does not touch
+ * component_wise() function of this class. Since it does not touch
* the order of indices within each, it may be worthwhile to first
* renumber using the Cuthill-McKee or a similar algorithm and
* afterwards renumbering component-wise. This will bring out the
* block.
*
*
- * @sect2{Cell-wise numbering for Discontinuous Galerkin FEM}
+ * <h3>Cell-wise numbering for Discontinuous Galerkin FEM</h3>
*
* One advantage of DGFEM is the fact, that it yields invertible
* blocks on the diagonal of the global matrix. these blocks are in
* preserved, so it may be useful to apply component_wise() first.
*
*
- * @sect2{Random renumbering}
+ * <h3>Random renumbering</h3>
*
- * The @p random function renumbers degrees of freedom randomly. This
+ * The random() function renumbers degrees of freedom randomly. This
* function is probably seldom of use, except to check the dependence of
* solvers (iterative or direct ones) on the numbering of the degrees
* of freedom. It uses the @p random_shuffle function from the C++
* standard library to do its work.
*
*
- * @sect2{Other renumberings}
- *
- * Apart from the ones discussed above, there are a number of other
- * renumbering schemes implemented in this class. Refer to the detailed
- * function listings for each of the functions of this class.
- *
- *
- * @sect2{Multigrid DoF numbering}
- *
- * Most algorithms also work on multigrid degree of freedom numberings. Refer
- * to the actual function declarations to get more information on this.
+ * <h3>Multigrid DoF numbering</h3>
*
+ * Most algorithms also work on multigrid degree of freedom
+ * numberings. Refer to the actual function declarations to get more
+ * information on this.
*
+ * @ingroup dofs
* @author Wolfgang Bangerth, Guido Kanschat, 1998, 1999, 2000, 2004
*/
class DoFRenumbering
{
public:
/**
- * Renumber the degrees of freedom
- * according to the Cuthill-McKee method,
- * eventually using the reverse numbering
- * scheme.
+ * Renumber the degrees of
+ * freedom according to the
+ * Cuthill-McKee method,
+ * eventually using the reverse
+ * numbering scheme.
*
- * See the general documentation of
- * this class for details on the
- * different methods.
+ * See the general documentation
+ * of this class for details on
+ * the different methods.
*/
template <int dim>
static void
/**
* Computes the renumbering
* vector needed by the
- * @p Cuthill_McKee
- * function. Does not perform the
- * renumbering on the
- * @p DoFHandler dofs but
+ * Cuthill_McKee() function. Does
+ * not perform the renumbering on
+ * the DoFHandler dofs but
* returns the renumbering
* vector.
*/
const std::vector<unsigned int> &starting_indices = std::vector<unsigned int>());
/**
- * Renumber the degrees of freedom
- * according to the Cuthill-McKee method,
- * eventually using the reverse numbering
- * scheme, in this case for a multigrid
- * numbering of degrees of freedom.
+ * Renumber the degrees of
+ * freedom according to the
+ * Cuthill-McKee method,
+ * eventually using the reverse
+ * numbering scheme, in this case
+ * for a multigrid numbering of
+ * degrees of freedom.
*
- * You can give a triangulation level to
- * which this function is to be applied.
- * Since with a level-wise numbering there
- * are no hanging nodes, no constraints
- * can be used, so the respective
- * parameter of the previous function is
+ * You can give a triangulation
+ * level to which this function
+ * is to be applied. Since with
+ * a level-wise numbering there
+ * are no hanging nodes, no
+ * constraints can be used, so
+ * the respective parameter of
+ * the previous function is
* ommitted.
*
- * See the general documentation of
- * this class for details on the
- * different methods.
+ * See the general documentation
+ * of this class for details on
+ * the different methods.
*/
template <int dim>
static void
/**
* Sort the degrees of freedom by
- * vector component. The numbering within
- * each component is not touched,
- * so a degree of freedom with index
- * $i$, belonging to some component,
- * and another degree of freedom
- * with index $j$ belonging to the same
+ * vector component. The
+ * numbering within each
+ * component is not touched, so a
+ * degree of freedom with index
+ * $i$, belonging to some
+ * component, and another degree
+ * of freedom with index $j$
+ * belonging to the same
* component will be assigned new
* indices $n(i)$ and $n(j)$ with
* $n(i)<n(j)$ if $i<j$ and
* You can specify that the
* components are ordered in a
* different way than suggested
- * by the @p FESystem object you
+ * by the FESystem object you
* use. To this end, set up the
* vector @p target_component
* such that the entry at index
- * @p i denotes the number of
- * the target component for dofs
- * with component @p i in the
- * @p FESystem. Naming the same
+ * @p i denotes the number of the
+ * target component for dofs with
+ * component @p i in the
+ * FESystem. Naming the same
* component more than once is
* possible and results in a
* blocking of several components
- * into one.
- * If you omit this argument, the
- * same order as given by the
- * finite element is used.
+ * into one. If you omit this
+ * argument, the same order as
+ * given by the finite element is
+ * used.
*
* If one of the base finite
* elements from which the global
/**
* Computes the renumbering
* vector needed by the
- * @p component_wise
- * functions. Does not perform the
- * renumbering on the
- * @p DoFHandler dofs but
- * returns the renumbering
- * vector.
+ * component_wise()
+ * functions. Does not perform
+ * the renumbering on the
+ * DoFHandler dofs but returns
+ * the renumbering vector.
*/
template <int dim, class ITERATOR, class ENDITERATOR>
static unsigned int
* Cell-wise renumbering for DG
* elements. This function takes
* the ordered set of cells in
- * @p cell_order, and makes sure
- * that all degrees of freedom in
- * a cell with higher index are
- * behind all degrees of freedom
- * of a cell with lower
- * index. The order inside a cell
- * bloock will be the same as
- * before this renumbering.
+ * <tt>cell_order</tt>, and makes
+ * sure that all degrees of
+ * freedom in a cell with higher
+ * index are behind all degrees
+ * of freedom of a cell with
+ * lower index. The order inside
+ * a cell bloock will be the same
+ * as before this renumbering.
*
* This function only works with
* Discontinuous Galerkin Finite
/**
* Computes the renumbering
* vector needed by the
- * @p cell_wise_dg
- * function. Does not perform the
- * renumbering on the
- * @p DoFHandler dofs but
+ * cell_wise_dg() function. Does
+ * not perform the renumbering on
+ * the DoFHandler dofs but
* returns the renumbering
* vector.
*/
* This function produces a
* downstream ordering of the
* mesh cells and calls
- * cell_wise_dg().
- * Therefore, it only works with
- * Discontinuous Galerkin Finite
- * Elements, i.e. all degrees of
- * freedom have to be associated
- * with the interior of the cell.
+ * cell_wise_dg(). Therefore, it
+ * only works with Discontinuous
+ * Galerkin Finite Elements,
+ * i.e. all degrees of freedom
+ * have to be associated with the
+ * interior of the cell.
*/
template <int dim>
static void
/**
* Computes the renumbering
* vector needed by the
- * @p downstream_dg
- * function. Does not perform the
- * renumbering on the
- * @p DoFHandler dofs but
+ * downstream_dg() function. Does
+ * not perform the renumbering on
+ * the DoFHandler dofs but
* returns the renumbering
* vector.
*/
* (counter)clockwise ordering of
* the mesh cells with respect to
* the hub @p center and calls
- * cell_wise_dg().
- * Therefore, it only works with
- * Discontinuous Galerkin Finite
- * Elements, i.e. all degrees of
- * freedom have to be associated
- * with the interior of the cell.
+ * cell_wise_dg(). Therefore, it
+ * only works with Discontinuous
+ * Galerkin Finite Elements,
+ * i.e. all degrees of freedom
+ * have to be associated with the
+ * interior of the cell.
*/
template <int dim>
static void
/**
* Computes the renumbering
* vector needed by the
- * @p clockwise_dg
- * functions. Does not perform the
- * renumbering on the
- * @p DoFHandler dofs but
+ * clockwise_dg() functions. Does
+ * not perform the renumbering on
+ * the DoFHandler dofs but
* returns the renumbering
* vector.
*/
/**
* Computes the renumbering
* vector needed by the
- * @p sort_selected_dofs_back
+ * sort_selected_dofs_back()
* function. Does not perform the
- * renumbering on the
- * @p DoFHandler dofs but
- * returns the renumbering
- * vector.
+ * renumbering on the DoFHandler
+ * dofs but returns the
+ * renumbering vector.
*/
template <int dim>
static void
/**
* Computes the renumbering
- * vector needed by the
- * @p random function. Does not
- * perform the renumbering on the
- * @p DoFHandler dofs but
- * returns the renumbering
- * vector.
+ * vector needed by the random()
+ * function. Does not perform the
+ * renumbering on the DoFHandler
+ * dofs but returns the
+ * renumbering vector.
*/
template <int dim>
static void
const DoFHandler<dim> &dof_handler);
/**
- * Renumber the degrees of freedom such
- * that they are associated with the
- * subdomain id of the cells they are
- * living on, i.e. first all degrees of
- * freedom that belong to cells with
- * subdomain zero, then all with
- * subdomain one, etc. This is useful
- * when doing parallel computations after
- * assigning subdomain ids using a
- * partitioner (see the
- * @p GridTools::partition_triangulation
+ * Renumber the degrees of
+ * freedom such that they are
+ * associated with the subdomain
+ * id of the cells they are
+ * living on, i.e. first all
+ * degrees of freedom that belong
+ * to cells with subdomain zero,
+ * then all with subdomain one,
+ * etc. This is useful when doing
+ * parallel computations after
+ * assigning subdomain ids using
+ * a partitioner (see the
+ * GridTools::partition_triangulation
* function for this).
*
* Note that degrees of freedom
- * associated with faces, edges, and
- * vertices may be associated with
- * multiple subdomains if they are
- * sitting on partition boundaries. It
- * would therefore be undefined with
- * which subdomain they have to be
- * associated. For this, we use what we
- * get from the
- * @p DoFTools::get_subdomain_association
+ * associated with faces, edges,
+ * and vertices may be associated
+ * with multiple subdomains if
+ * they are sitting on partition
+ * boundaries. It would therefore
+ * be undefined with which
+ * subdomain they have to be
+ * associated. For this, we use
+ * what we get from the
+ * DoFTools::get_subdomain_association
* function.
*
* The algorithm is stable, i.e. if
subdomain_wise (DoFHandler<dim> &dof_handler);
/**
- * Computes the renumbering vector needed
- * by the @p subdomain_wise
+ * Computes the renumbering
+ * vector needed by the
+ * subdomain_wise()
* function. Does not perform the
- * renumbering on the @p DoFHandler dofs
- * but returns the renumbering vector.
+ * renumbering on the @p
+ * DoFHandler dofs but returns
+ * the renumbering vector.
*/
template <int dim>
static void
/**
* Exception
+ *
+ * @ingroup Exceptions
*/
DeclException0 (ExcRenumberingIncomplete);
/**
* Exception
+ *
+ * @ingroup Exceptions
*/
DeclException0 (ExcInvalidComponentOrder);
/**
- * Exception. The function is
- * only implemented for
- * Discontinuous Galerkin Finite
- * elements.
+ * The function is only
+ * implemented for Discontinuous
+ * Galerkin Finite elements.
+ *
+ * @ingroup Exceptions
*/
DeclException0 (ExcNotDGFEM);
};
* object of class DoFTools.
*
*
- * @sect3{Setting up sparsity patterns}
+ * <h3>Setting up sparsity patterns</h3>
*
* When assembling system matrices, the entries are usually of the form
* $a_{ij} = a(\phi_i, \phi_j)$, where $a$ is a bilinear functional, often an
*
*
*
- * @sect3{DoF numberings on boundaries}
+ * <h3>DoF numberings on boundaries</h3>
*
* When projecting the traces of functions to the boundary or parts
* thereof, one needs to build matrices and vectors with the degrees
* of freedom on the boundary. What is needed in this case is a
* numbering of the boundary degrees of freedom, starting from zero on
* and not considering the degrees of freedom in the interior. The
- * @p map_dof_to_boundary_indices function does exactly this, by
+ * map_dof_to_boundary_indices() function does exactly this, by
* providing a vector with as many entries as there are degrees of
* freedom on the whole domain, with each entry being the number in
* the numbering of the boundary or
- * DoFHandler@p ::invalid_dof_index if the dof is not on the
+ * DoFHandler::invalid_dof_index if the dof is not on the
* boundary. You should always use this function to get the mapping
* between local (boundary) and the global numbers, for example to
* build the mass matrix on the boundary, or to get the global index
* algorithm, you are better off if you just accept the mapping `as
* is'.
*
- * Actually, there are two @p map_dof_to_boundary_indices functions,
+ * Actually, there are two map_dof_to_boundary_indices() functions,
* one producing a numbering for all boundary degrees of freedom and
* one producing a numbering for only parts of the boundary, namely
* those parts for which the boundary indicator is listed in a set of
* indices of degrees of freedom on different parts may be intermixed.
*
* Degrees of freedom on the boundary but not on one of the specified
- * boundary parts are given the index @p invalid_dof_index, as if
+ * boundary parts are given the index #invalid_dof_index, as if
* they were in the interior. If no boundary indicator was given or if
* no face of a cell has a boundary indicator contained in the given
* list, the vector of new indices consists solely of
- * @p invalid_dof_indexs.
+ * #invalid_dof_index.
*
* The question what a degree of freedom on the boundary is, is not so
* easy. It should really be a degree of freedom of which the
* least for Lagrange elements this definition is equal to the
* statement that the off-point of the shape function, i.e. the point
* where the function assumes its nominal value (for Lagrange elements
- * this is the point where it has the function value @p 1), is
+ * this is the point where it has the function value 1), is
* located on the boundary. We do not check this directly, the
* criterion is rather defined through the information the finite
* element class gives: the FiniteElementBase class defines the
* is, is a secret of the finite element (well, you can ask it, but we
* don't do it here) and not relevant in this context.
*
- *
+ * @ingroup dofs
* @author Wolfgang Bangerth, Guido Kanschat and others, 1998 - 2005
*/
class DoFTools
*/
nonzero
};
+ /**
+ * @name Sparsity Pattern Generation
+ * @{
+ */
+ /**
+ * On initializing a
+ * SparsityPattern, the number of
+ * entries required for each row
+ * must be known. This function
+ * estimates the maximum number
+ * of coupling degrees of freedom
+ * for each row of the sparsity
+ * pattern.
+ *
+ * In systems of equations,
+ * couplings between components
+ * <b>not</b> coupling in the
+ * differential equation can be
+ * eliminated by the two optional
+ * tables.
+ *
+ * @param dofs The DoFHandler
+ * @param row_lengths The vector
+ * containing the resulting row
+ * lengths. It must have the
+ * length DoFHandler::n_dofs().
+ *
+ * @param couplings
+ * Optional argument for
+ * optimizing out non-coupled
+ * components of a system. If the
+ * entry at position <i>(i,j)</i>
+ * is #none, then couplings
+ * between the related dofs are
+ * neglected. See Coupling.
+ *
+ * @param flux_couplings:
+ * similar to the previous
+ * argument allows neglecting
+ * couplings. Here, couplings due
+ * to flux operators on faces are
+ * considered. See #Coupling.
+ */
+ template<int dim>
+ static
+ void compute_row_length_vector(
+ const DoFHandler<dim>& dofs,
+ std::vector<unsigned int>& row_lengths,
+ Table<2,Coupling> couplings /*= typename Table<2,Coupling>()*/,
+ Table<2,Coupling> flux_couplings /*= Table<2,Coupling>()*/);
/**
* Locate non-zero entries of the
SparsityPattern &sparsity,
const Table<2,Coupling>& int_mask,
const Table<2,Coupling>& flux_mask);
+
+ //@}
+ /**
+ * @name Hanging Nodes
+ * @{
+ */
/**
* Make up the constraints which
static void
make_hanging_node_constraints (const DoFHandler<3> &dof_handler,
ConstraintMatrix &constraints);
+ //@}
/**
* Take a vector of values which live on
* note that the resulting field will not
* be continuous at hanging nodes. This
* can, however, easily be arranged by
- * calling the appropraite @p distribute
+ * calling the appropriate @p distribute
* function of a ConstraintMatrix
* object created for this
* DoFHandler object, after the
* It is assumed that the output
* vector @p dof_data already
* has the right size,
- * i.e. <tt>n_dofs()</tt> elements.
+ * i.e. n_dofs() elements.
*
* This function cannot be used
* if the finite element in use
const std::vector<bool> &component_select,
std::vector<bool> &selected_dofs,
const std::set<unsigned char> &boundary_indicators = std::set<unsigned char>());
-
+ /**
+ * @name Hanging Nodes
+ * @{
+ */
+
/**
* Select all dofs that will be
* constrained by interface
static void
extract_hanging_node_dofs (const DoFHandler<3> &dof_handler,
std::vector<bool> &selected_dofs);
-
+ //@}
+
/**
* Flag all those degrees of
* freedom which are on cells
/**
* Declare the type as discussed
* above. Since we can't name it
- * @p FunctionMap (as that would
+ * FunctionMap (as that would
* ambiguate a possible
* constructor of this class),
* name it in the fashion of the
* continuously differentiable
* over cell boundaries.
*
- * The value <i>L<sup>2</sup></i>
+ * <li> <i>L<sup>2</sup></i>
* indicates that the element is
* discontinuous. Since
* discontinuous elements have no
* special way in the sense that
* it is <b>not</b> implied by
* any higher conformity.
+ * </ol>
*
* In order to test if a finite
* element conforms to a certain
* @param degree
* Maximal polynomial degree in a
* single direction.
+ * @param conformity The finite
+ * element space has continuity
+ * of this Sobolev space.
*/
FiniteElementData (const std::vector<unsigned int> &dofs_per_object,
const unsigned int n_components,
#include <numeric>
+//TODO:[GK] Traverse faces only once using flags
+
+//TODO:[GK] This function is not finished yet!!!
+template <int dim>
+void
+DoFTools::compute_row_length_vector(
+ const DoFHandler<dim>& dofs,
+ std::vector<unsigned int>& row_lengths,
+ Table<2,Coupling> couplings,
+ Table<2,Coupling> flux_couplings)
+{
+ const FiniteElement<dim>& fe = dofs.get_fe();
+ const unsigned int ncomp = fe.n_components;
+
+ Assert (row_lengths.size() == dofs.n_dofs(),
+ ExcDimensionMismatch(row_lengths.size(), dofs.n_dofs()));
+
+ Assert (couplings.n_rows() == couplings.n_cols(),
+ ExcDimensionMismatch(couplings.n_rows(), couplings.n_cols()));
+ Assert (flux_couplings.n_rows() == flux_couplings.n_cols(),
+ ExcDimensionMismatch(flux_couplings.n_rows(), flux_couplings.n_cols()));
+
+ if (couplings.n_rows() > 1)
+ Assert (couplings.n_rows() == ncomp,
+ ExcDimensionMismatch(couplings.n_rows(), ncomp));
+ if (flux_couplings.n_rows() > 1)
+ Assert (flux_couplings.n_rows() == ncomp,
+ ExcDimensionMismatch(flux_couplings.n_rows(), ncomp));
+
+ // Function starts here by
+ // resetting the counters.
+ std::fill(row_lengths.begin(), row_lengths.end(), 0);
+
+ const typename DoFHandler<dim>::cell_iterator end = dofs.end();
+ typename DoFHandler<dim>::active_cell_iterator cell;
+ std::vector<unsigned int> indices(fe.dofs_per_cell);
+
+ for (cell = dofs.begin(); cell != end; ++cell)
+ {
+ indices.resize(fe.dofs_per_cell);
+ cell->get_dof_indices(indices);
+ // First the simple case where
+ // all degrees of freedom on a
+ // cell couple.
+ if (couplings.n_rows() <= 1)
+ {
+ unsigned int i = 0;
+ // First, dofs on
+ // vertices. We assume that
+ // each vertex dof couples
+ // with all dofs on
+ // adjacent grid cells.
+
+ // Since two cells share a
+ // facein 2d, we can subtract
+ // the dofs per face,
+ // unless we have hanging
+ // nodes.
+
+ // In 3d, 8 cells share 12
+ // faces.
+ unsigned int increment = fe.dofs_per_cell - fe.dofs_per_face;
+ while (i < fe.first_line_index)
+ row_lengths[indices[i++]] += increment;
+ // Lines are already cells
+ // in 1D
+ increment = (dim>1)
+ ? fe.dofs_per_cell - fe.dofs_per_face
+ : fe.dofs_per_cell;
+ while (i < fe.first_quad_index)
+ row_lengths[indices[i++]] += increment;
+ // Now quads in 2D and 3D
+ increment = (dim>2)
+ ? fe.dofs_per_cell - fe.dofs_per_face
+ : fe.dofs_per_cell;
+ while (i < fe.first_hex_index)
+ row_lengths[indices[i++]] += increment;
+ // Finally, cells in 3D
+ increment = fe.dofs_per_cell;
+ while (i < fe.ofs_per_cell)
+ row_lengths[indices[i++]] += increment;
+ }
+ // At this point, we have
+ // counted all dofs
+ // contributiong from cells
+ // coupled topologically to the
+ // adjacent cells, but we
+ // subtracted some faces.
+
+ // Now, let's go by the faces
+ // and add the missing
+ // contribution as well as the
+ // flux contributions.
+ for (unsigned int iface=0;iface<GeometryInfo<dim>::faces_per_cell;++iface)
+ {
+ typename DoFHandler<dim>::face_iterator face = cell->face(iface);
+
+ indices.resize(fe.dofs_per_face);
+ face->get_dof_indices(indices);
+ for (i=0;i<fe.dofs_per_face;++i)
+ row_lengths[indices[i]] += fe.dofs_per_face;
+ }
+ }
+}
+
+
template <int dim, class SparsityPattern>
void
DoFTools::make_sparsity_pattern (
projects / dealii-svn.git / commitdiff