From c2c7912425b4960c3ca701f8813b471b42d56a0c Mon Sep 17 00:00:00 2001 From: guido Date: Thu, 30 Jun 2005 21:23:46 +0000 Subject: [PATCH] new function DoFTools::compute_row_length_vector Doxygen group dofs filled git-svn-id: https://svn.dealii.org/trunk@11023 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/dofs/dof_accessor.h | 76 +-- .../deal.II/include/dofs/dof_constraints.h | 27 +- deal.II/deal.II/include/dofs/dof_handler.h | 118 ++++- deal.II/deal.II/include/dofs/dof_levels.h | 6 +- .../deal.II/include/dofs/dof_renumbering.h | 447 +++++++++--------- deal.II/deal.II/include/dofs/dof_tools.h | 88 +++- deal.II/deal.II/include/dofs/function_map.h | 2 +- deal.II/deal.II/include/fe/fe_base.h | 6 +- deal.II/deal.II/source/dofs/dof_tools.cc | 106 +++++ 9 files changed, 598 insertions(+), 278 deletions(-) diff --git a/deal.II/deal.II/include/dofs/dof_accessor.h b/deal.II/deal.II/include/dofs/dof_accessor.h index b5b5b61d7d..903ac25a6b 100644 --- a/deal.II/deal.II/include/dofs/dof_accessor.h +++ b/deal.II/deal.II/include/dofs/dof_accessor.h @@ -14,9 +14,6 @@ #define __deal2__dof_accessor_h -/*------------------------- dof_iterator.h ------------------------*/ - - #include #include #include @@ -46,8 +43,8 @@ template class DoFObjectAccessor<3, dim>; /** * 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 @@ -63,8 +60,9 @@ template class DoFObjectAccessor<3, dim>; * 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 @@ -78,7 +76,7 @@ class DoFAccessor /** * 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 @@ -117,22 +115,34 @@ class DoFAccessor /** * 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); @@ -158,6 +168,7 @@ class DoFAccessor * 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 @@ -184,6 +195,7 @@ class DoFObjectAccessor_Inheritance * 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 @@ -219,7 +231,7 @@ class DoFObjectAccessor_Inheritance * lines in 1D-, 2D-, etc dimensions). * * - * @sect3{Usage} + *

Usage

* * The DoFDimensionInfo classes inherited by the DoFHandler classes * declare typedefs to iterators using the accessors declared in this class @@ -228,12 +240,12 @@ class DoFObjectAccessor_Inheritance * as they provide easier typing (much less complicated names!). * * - * @sect3{Notes about the class hierarchy structure} + *

Notes about the class hierarchy structure

* * 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 @@ -272,27 +284,27 @@ class DoFObjectAccessor : public DoFAccessor, const AccessorData *local_data); /** - * Index of the @p ith degree + * Index of the ith 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 ith 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 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 degree * on the @p vertexth vertex to @p index. */ void set_vertex_dof_index (const unsigned int vertex, @@ -343,7 +355,7 @@ class DoFObjectAccessor : public DoFAccessor, /** * 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 @@ -380,21 +392,21 @@ class DoFObjectAccessor : public DoFAccessor, OutputVector &values) const; /** - * Pointer to the @p ith line + * Pointer to the ith line * bounding this Object. */ TriaIterator > line (const unsigned int i) const; /** - * Pointer to the @p ith quad + * Pointer to the ith quad * bounding this Object. */ TriaIterator > quad (const unsigned int i) const; /** - * @p ith child as a DoFObjectAccessor + * ith child as a DoFObjectAccessor * iterator. This function is needed since * the child function of the base * class returns a hex accessor without @@ -479,7 +491,7 @@ class DoFObjectAccessor<0, dim> : public DoFAccessor, * lines in 1D-, 2D-, etc dimensions). * * - * @sect3{Usage} + *

Usage

* * The DoFDimensionInfo classes inherited by the DoFHandler classes * declare typedefs to iterators using the accessors declared in this class @@ -488,7 +500,7 @@ class DoFObjectAccessor<0, dim> : public DoFAccessor, * as they provide easier typing (much less complicated names!). * * - * @sect3{Notes about the class hierarchy structure} + *

Notes about the class hierarchy structure

* * Inheritance from DoFObjectAccessor_Inheritance<1,dim>::BaseClass yields * inheritance from CellAccessor<1> if dim==1 and from @@ -496,6 +508,7 @@ class DoFObjectAccessor<0, dim> : public DoFAccessor, * 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 @@ -601,7 +614,7 @@ class DoFObjectAccessor<1, dim> : public DoFAccessor, /** * 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 @@ -689,7 +702,7 @@ class DoFObjectAccessor<1, dim> : public DoFAccessor, /** * Grant access to the degrees of freedom located on quads. * - * @ref DoFObjectAccessor + * @ingroup dofs */ template class DoFObjectAccessor<2, dim> : public DoFAccessor, @@ -795,7 +808,7 @@ class DoFObjectAccessor<2, dim> : public DoFAccessor, /** * 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 @@ -891,7 +904,7 @@ class DoFObjectAccessor<2, dim> : public DoFAccessor, /** * Grant access to the degrees of freedom located on hexes. * - * @ref DoFObjectAccessor + * @ingroup dofs */ template class DoFObjectAccessor<3, dim> : public DoFAccessor, @@ -997,7 +1010,7 @@ class DoFObjectAccessor<3, dim> : public DoFAccessor, /** * 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 @@ -1115,6 +1128,7 @@ class DoFObjectAccessor<3, dim> : public DoFAccessor, * CellAccessor, which makes the functions of this class available to the * DoFCellAccessor class as well. * + * @ingroup dofs * @author Wolfgang Bangerth, 1998 */ template @@ -1237,7 +1251,7 @@ class DoFCellAccessor : public DoFObjectAccessor * 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 @@ -1253,7 +1267,7 @@ class DoFCellAccessor : public DoFObjectAccessor /** * 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 @@ -1265,7 +1279,7 @@ class DoFCellAccessor : public DoFObjectAccessor * 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 @@ -1273,7 +1287,7 @@ class DoFCellAccessor : public DoFObjectAccessor * 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 @@ -1319,7 +1333,7 @@ class DoFCellAccessor : public DoFObjectAccessor * 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 diff --git a/deal.II/deal.II/include/dofs/dof_constraints.h b/deal.II/deal.II/include/dofs/dof_constraints.h index 1fe8dd24e2..8597932fce 100644 --- a/deal.II/deal.II/include/dofs/dof_constraints.h +++ b/deal.II/deal.II/include/dofs/dof_constraints.h @@ -75,11 +75,11 @@ class BlockIndices; * which has to applied to the solution vector. * * - * @sect3{Condensing matrices and sparsity patterns} + *

Condensing matrices and sparsity patterns

* * 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 @@ -147,7 +147,7 @@ class BlockIndices; * in the already condensed form. * * - * @sect3{Condensing vectors} + *

Condensing vectors

* * 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 @@ -161,7 +161,7 @@ class BlockIndices; * techniques as mentioned above to avoid their use. * * - * @sect3{Avoiding explicit condensation} + *

Avoiding explicit condensation

* * Sometimes, one wants to avoid condensation at all. This may be the case * since condensation is an expensive operation, or because no condense() @@ -181,7 +181,7 @@ class BlockIndices; * system. * * - * @sect3{Distributing constraints} + *

Distributing constraints

* * After solving the condensed system of equations, the solution vector has to * be redistributed. This is done by the two @p distribute function, one @@ -194,6 +194,7 @@ class BlockIndices; * 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 @@ -793,14 +794,20 @@ class ConstraintMatrix : public Subscriptor /** * Exception + * + * @ingroup Exceptions */ DeclException0 (ExcMatrixIsClosed); /** * Exception + * + * @ingroup Exceptions */ DeclException0 (ExcMatrixNotClosed); /** * Exception + * + * @ingroup Exceptions */ DeclException1 (ExcLineInexistant, unsigned int, @@ -808,10 +815,14 @@ class ConstraintMatrix : public Subscriptor << " does not exist."); /** * Exception + * + * @ingroup Exceptions */ DeclException0 (ExcWrongDimension); /** * Exception + * + * @ingroup Exceptions */ DeclException4 (ExcEntryAlreadyExists, int, int, double, double, @@ -821,6 +832,8 @@ class ConstraintMatrix : public Subscriptor << "by " << (arg4-arg3) << "."); /** * Exception + * + * @ingroup Exceptions */ DeclException2 (ExcDoFConstrainedToConstrainedDoF, int, int, @@ -829,6 +842,8 @@ class ConstraintMatrix : public Subscriptor << ", but that one is also constrained. This is not allowed!"); /** * Exception. + * + * @ingroup Exceptions */ DeclException1 (ExcDoFIsConstrainedFromBothObjects, int, @@ -836,6 +851,8 @@ class ConstraintMatrix : public Subscriptor << " is constrained from both object in a merge operation."); /** * Exception + * + * @ingroup Exceptions */ DeclException1 (ExcDoFIsConstrainedToConstrainedDoF, int, diff --git a/deal.II/deal.II/include/dofs/dof_handler.h b/deal.II/deal.II/include/dofs/dof_handler.h index b488fe7f58..30d68b5ba5 100644 --- a/deal.II/deal.II/include/dofs/dof_handler.h +++ b/deal.II/deal.II/include/dofs/dof_handler.h @@ -45,11 +45,111 @@ template class Triangulation; * @ref DoFDimensionInfo<1> * @ref DoFDimensionInfo<2> * + * @ingroup dofs * @author Wolfgang Bangerth, 1998 */ template 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 + *
    + *
  1. #line_iterator + *
  2. #quad_iterator + *
  3. #hex_iterator + *
+ */ + typedef void * cell_iterator; + /** + * All active cells of the DoFHandler + * Depending on the dimension + * this is typedefed to + *
    + *
  1. #active_line_iterator + *
  2. #active_quad_iterator + *
  3. #active_hex_iterator + *
+ */ + 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 + *
    + *
  1. void* + *
  2. #line_iterator + *
  3. #quad_iterator + *
+ */ + typedef void * face_iterator; + /** + * All active faces of the DoFHandler + * Depending on the dimension + * this is typedefed to + *
    + *
  1. void* + *
  2. #active_line_iterator + *
  3. #active_quad_iterator + *
+ */ + typedef void * active_face_iterator; }; @@ -59,6 +159,7 @@ class DoFDimensionInfo * * The types have the same meaning as those declared in TriaDimensionInfo<2>. * + * @ingroup dofs * @author Wolfgang Bangerth, 1998 */ template <> @@ -93,6 +194,7 @@ class DoFDimensionInfo<1> * * The types have the same meaning as those declared in TriaDimensionInfo<2>. * + * @ingroup dofs * @author Wolfgang Bangerth, 1998 */ template <> @@ -127,6 +229,7 @@ class DoFDimensionInfo<2> * * The types have the same meaning as those declared in TriaDimensionInfo<3>. * + * @ingroup dofs * @author Wolfgang Bangerth, 1998 */ template <> @@ -184,10 +287,10 @@ class DoFDimensionInfo<3> * global matrices as well as matrices living on (parts of) the boundary. * * - * @sect3{Distribution of indices for degrees of freedom} + *

Distribution of indices for degrees of freedom

* * The degrees of freedom (`dofs') are distributed on the given triangulation - * by the function distribute_dofs(). 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 @@ -206,13 +309,13 @@ class DoFDimensionInfo<3> * the implemented algorithms. * * - * @sect3{User defined renumbering schemes} + *

User defined renumbering schemes

* * 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 - * renumber_dofs(vector) 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 @@ -223,8 +326,7 @@ class DoFDimensionInfo<3> * renumber_dofs(vector) with the array, which converts old * into new degree of freedom indices. * - * - * + * @ingroup dofs * @author Wolfgang Bangerth, 1998 */ template @@ -1170,12 +1272,12 @@ class DoFHandler : public Subscriptor, */ std::vector vertex_dofs; - /** + /* * Make accessor objects friends. */ template friend class DoFAccessor; - /** + /* * Make accessor objects friends. */ template friend class DoFObjectAccessor; diff --git a/deal.II/deal.II/include/dofs/dof_levels.h b/deal.II/deal.II/include/dofs/dof_levels.h index bf4bf4a047..0b85247e2f 100644 --- a/deal.II/deal.II/include/dofs/dof_levels.h +++ b/deal.II/deal.II/include/dofs/dof_levels.h @@ -24,6 +24,7 @@ * class, but do not actually use it. Rather, only specializations of * this class are used. * + * @ingroup dofs * @author Wolfgang Bangerth, 1998 */ template @@ -43,7 +44,7 @@ class DoFLevel * Store the indices of the degrees of freedom which are located on * the lines. * - * @sect3{Information for all DoFLevel classes} + *

Information for all DoFLevel classes

* * The DoFLevel classes * store the global indices of the degrees of freedom for each cell on a @@ -86,6 +87,7 @@ class DoFLevel * is used, the indices are stored in the @p vertex_dofs array of the * DoFHandler class. * + * @ingroup dofs * @author Wolfgang Bangerth, 1998 */ template <> @@ -114,6 +116,7 @@ class DoFLevel<1> * 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 <> @@ -142,6 +145,7 @@ class DoFLevel<2> : public DoFLevel<1> * 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 <> diff --git a/deal.II/deal.II/include/dofs/dof_renumbering.h b/deal.II/deal.II/include/dofs/dof_renumbering.h index 7b5f601ab9..2a8aa51a85 100644 --- a/deal.II/deal.II/include/dofs/dof_renumbering.h +++ b/deal.II/deal.II/include/dofs/dof_renumbering.h @@ -27,100 +27,107 @@ * 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 + *

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 * 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. + *

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::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. + *

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. * * 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 @@ -145,7 +152,7 @@ * may be difficult, however, and in many cases will not justify the effort. * * - * @sect2{Component-wise numbering} + *

Component-wise numbering

* * For finite elements composed of several base elements using the FESystem * class, or for elements which provide several components themselves, it @@ -155,7 +162,7 @@ * 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 @@ -163,7 +170,7 @@ * block. * * - * @sect2{Cell-wise numbering for Discontinuous Galerkin FEM} + *

Cell-wise numbering for Discontinuous Galerkin FEM

* * One advantage of DGFEM is the fact, that it yields invertible * blocks on the diagonal of the global matrix. these blocks are in @@ -179,42 +186,37 @@ * preserved, so it may be useful to apply component_wise() first. * * - * @sect2{Random renumbering} + *

Random renumbering

* - * 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. + *

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. * + * @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 static void @@ -226,10 +228,9 @@ class DoFRenumbering /** * 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. */ @@ -242,23 +243,27 @@ class DoFRenumbering const std::vector &starting_indices = std::vector()); /** - * 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 static void @@ -269,12 +274,14 @@ class DoFRenumbering /** * 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) static unsigned int @@ -361,14 +367,14 @@ class DoFRenumbering * 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. + * 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. * * This function only works with * Discontinuous Galerkin Finite @@ -384,10 +390,9 @@ class DoFRenumbering /** * 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. */ @@ -430,12 +435,12 @@ class DoFRenumbering * 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 static void @@ -458,10 +463,9 @@ class DoFRenumbering /** * 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. */ @@ -479,12 +483,12 @@ class DoFRenumbering * (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 static void @@ -507,10 +511,9 @@ class DoFRenumbering /** * 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. */ @@ -542,12 +545,11 @@ class DoFRenumbering /** * 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 static void @@ -565,12 +567,11 @@ class DoFRenumbering /** * 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 static void @@ -578,29 +579,32 @@ class DoFRenumbering const DoFHandler &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 @@ -614,11 +618,13 @@ class DoFRenumbering subdomain_wise (DoFHandler &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 static void @@ -627,17 +633,22 @@ class DoFRenumbering /** * 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); }; diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index 0145f36989..10a53886b7 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -46,7 +46,7 @@ template class Mapping; * object of class DoFTools. * * - * @sect3{Setting up sparsity patterns} + *

Setting up sparsity patterns

* * 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 @@ -70,18 +70,18 @@ template class Mapping; * * * - * @sect3{DoF numberings on boundaries} + *

DoF numberings on boundaries

* * 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 @@ -99,7 +99,7 @@ template class Mapping; * 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 @@ -114,11 +114,11 @@ template class Mapping; * 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 @@ -126,7 +126,7 @@ template class Mapping; * 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 @@ -141,7 +141,7 @@ template class Mapping; * 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 @@ -173,7 +173,57 @@ 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 + * not 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,j) + * 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 + static + void compute_row_length_vector( + const DoFHandler& dofs, + std::vector& row_lengths, + Table<2,Coupling> couplings /*= typename Table<2,Coupling>()*/, + Table<2,Coupling> flux_couplings /*= Table<2,Coupling>()*/); /** * Locate non-zero entries of the @@ -502,6 +552,12 @@ class DoFTools SparsityPattern &sparsity, const Table<2,Coupling>& int_mask, const Table<2,Coupling>& flux_mask); + + //@} + /** + * @name Hanging Nodes + * @{ + */ /** * Make up the constraints which @@ -550,6 +606,7 @@ class DoFTools static void make_hanging_node_constraints (const DoFHandler<3> &dof_handler, ConstraintMatrix &constraints); + //@} /** * Take a vector of values which live on @@ -561,7 +618,7 @@ class DoFTools * 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 @@ -600,7 +657,7 @@ class DoFTools * It is assumed that the output * vector @p dof_data already * has the right size, - * i.e. n_dofs() elements. + * i.e. n_dofs() elements. * * This function cannot be used * if the finite element in use @@ -746,7 +803,11 @@ class DoFTools const std::vector &component_select, std::vector &selected_dofs, const std::set &boundary_indicators = std::set()); - + /** + * @name Hanging Nodes + * @{ + */ + /** * Select all dofs that will be * constrained by interface @@ -778,7 +839,8 @@ class DoFTools static void extract_hanging_node_dofs (const DoFHandler<3> &dof_handler, std::vector &selected_dofs); - + //@} + /** * Flag all those degrees of * freedom which are on cells diff --git a/deal.II/deal.II/include/dofs/function_map.h b/deal.II/deal.II/include/dofs/function_map.h index c72db9069d..0f42724ff9 100644 --- a/deal.II/deal.II/include/dofs/function_map.h +++ b/deal.II/deal.II/include/dofs/function_map.h @@ -39,7 +39,7 @@ struct FunctionMap /** * 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 diff --git a/deal.II/deal.II/include/fe/fe_base.h b/deal.II/deal.II/include/fe/fe_base.h index 0e7e8d9d27..e755bc5547 100644 --- a/deal.II/deal.II/include/fe/fe_base.h +++ b/deal.II/deal.II/include/fe/fe_base.h @@ -84,7 +84,7 @@ class FiniteElementData * continuously differentiable * over cell boundaries. * - * The value L2 + *
  • L2 * indicates that the element is * discontinuous. Since * discontinuous elements have no @@ -96,6 +96,7 @@ class FiniteElementData * special way in the sense that * it is not implied by * any higher conformity. + * * * In order to test if a finite * element conforms to a certain @@ -294,6 +295,9 @@ class FiniteElementData * @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 &dofs_per_object, const unsigned int n_components, diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index 085031b092..0111f6b00e 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -37,6 +37,112 @@ #include +//TODO:[GK] Traverse faces only once using flags + +//TODO:[GK] This function is not finished yet!!! +template +void +DoFTools::compute_row_length_vector( + const DoFHandler& dofs, + std::vector& row_lengths, + Table<2,Coupling> couplings, + Table<2,Coupling> flux_couplings) +{ + const FiniteElement& 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::cell_iterator end = dofs.end(); + typename DoFHandler::active_cell_iterator cell; + std::vector 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::faces_per_cell;++iface) + { + typename DoFHandler::face_iterator face = cell->face(iface); + + indices.resize(fe.dofs_per_face); + face->get_dof_indices(indices); + for (i=0;i void DoFTools::make_sparsity_pattern ( -- 2.39.5