From 2763d089986e4181308ca58e2e88279f2b0988b3 Mon Sep 17 00:00:00 2001 From: bangerth Date: Fri, 7 Mar 2008 22:11:09 +0000 Subject: [PATCH] Convert DoFRenumbering from a class with all static members to a namespace. git-svn-id: https://svn.dealii.org/trunk@15878 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/include/dofs/dof_renumbering.h | 1429 ++++--- .../deal.II/source/dofs/dof_renumbering.cc | 3498 +++++++++-------- deal.II/doc/news/changes.h | 18 +- 3 files changed, 2500 insertions(+), 2445 deletions(-) diff --git a/deal.II/deal.II/include/dofs/dof_renumbering.h b/deal.II/deal.II/include/dofs/dof_renumbering.h index 632ca14b95..a85a2bd144 100644 --- a/deal.II/deal.II/include/dofs/dof_renumbering.h +++ b/deal.II/deal.II/include/dofs/dof_renumbering.h @@ -367,169 +367,33 @@ DEAL_II_NAMESPACE_OPEN * @ingroup dofs * @author Wolfgang Bangerth, Guido Kanschat, 1998, 1999, 2000, 2004, 2007, 2008 */ -class DoFRenumbering +namespace DoFRenumbering { - public: - /** - * A namespace for the - * implementation of some - * renumbering algorithms based - * on algorithms implemented in - * the Boost Graph Library (BGL) - * by Jeremy Siek and others. - * - * While often slighty slower to - * compute, the algorithms using - * BOOST often lead to matrices - * with smaller bandwidths and - * sparse ILUs based on this - * numbering are therefore more - * efficient. - * - * For a comparison of these - * algorithms with the ones - * defined in DoFRenumbering, see - * the comparison section in the - * documentation of the - * DoFRenumbering namespace. - */ - class boost - { - public: - /** - * Renumber the degrees of - * freedom according to the - * Cuthill-McKee method, - * eventually using the reverse - * numbering scheme. - * - * See the general - * documentation of the - * parent class for details - * on the different methods. - * - * As an example of the - * results of this algorithm, - * take a look at the - * comparison of various - * algorithms in the - * documentation of the - * DoFRenumbering namespace. - */ - template - static void - Cuthill_McKee (DH& dof_handler, - const bool reversed_numbering = false, - const bool use_constraints = false); - - /** - * Computes the renumbering - * vector needed by the - * Cuthill_McKee() function. Does - * not perform the renumbering on - * the DoFHandler dofs but - * returns the renumbering - * vector. - */ - template - static void - compute_Cuthill_McKee (std::vector& new_dof_indices, - const DH&, - const bool reversed_numbering = false, - const bool use_constraints = false); - - /** - * Renumber the degrees of - * freedom based on the BOOST - * implementation of the King - * algorithm. This often - * results in slightly larger - * (by a few percent) - * bandwidths than the - * Cuthill-McKee algorithm, - * but sparse ILUs are often - * slightly (also by a few - * percent) better - * preconditioners. - * - * As an example of the - * results of this algorithm, - * take a look at the - * comparison of various - * algorithms in the - * documentation of the - * DoFRenumbering namespace. - * - * This algorithm is used in - * @ref step_22 "step-22". - */ - template - static void - king_ordering (DH& dof_handler, - const bool reversed_numbering = false, - const bool use_constraints = false); - - /** - * Compute the renumbering - * for the King algorithm but - * do not actually renumber - * the degrees of freedom in - * the DoF handler argument. - */ - template - static void - compute_king_ordering (std::vector& new_dof_indices, - const DH&, - const bool reversed_numbering = false, - const bool use_constraints = false); - - /** - * Renumber the degrees of - * freedom based on the BOOST - * implementation of the - * minimum degree - * algorithm. Unlike the - * Cuthill-McKee algorithm, - * this algorithm does not - * attempt to minimize the - * bandwidth of a matrix but - * to minimize the amount of - * fill-in when doing an LU - * decomposition. It may - * sometimes yield better - * ILUs because of this - * property. - * - * As an example of the - * results of this algorithm, - * take a look at the - * comparison of various - * algorithms in the - * documentation of the - * DoFRenumbering namespace. - */ - template - static void - minimum_degree (DH& dof_handler, - const bool reversed_numbering = false, - const bool use_constraints = false); - - /** - * Compute the renumbering - * for the minimum degree - * algorithm but do not - * actually renumber the - * degrees of freedom in the - * DoF handler argument. - */ - template - static void - compute_minimum_degree (std::vector& new_dof_indices, - const DH&, - const bool reversed_numbering = false, - const bool use_constraints = false); - }; - + /** + * A namespace for the + * implementation of some + * renumbering algorithms based + * on algorithms implemented in + * the Boost Graph Library (BGL) + * by Jeremy Siek and others. + * + * While often slighty slower to + * compute, the algorithms using + * BOOST often lead to matrices + * with smaller bandwidths and + * sparse ILUs based on this + * numbering are therefore more + * efficient. + * + * For a comparison of these + * algorithms with the ones + * defined in DoFRenumbering, see + * the comparison section in the + * documentation of the + * DoFRenumbering namespace. + */ + namespace boost + { /** * Renumber the degrees of * freedom according to the @@ -537,23 +401,24 @@ class DoFRenumbering * eventually using the reverse * numbering scheme. * - * See the general documentation - * of this class for details on - * the different methods. + * See the general + * documentation of the + * parent class for details + * on the different methods. * - * As an example of the results - * of this algorithm, take a look - * at the comparison of various + * As an example of the + * results of this algorithm, + * take a look at the + * comparison of various * algorithms in the * documentation of the * DoFRenumbering namespace. */ template - static void + void Cuthill_McKee (DH& dof_handler, const bool reversed_numbering = false, - const bool use_constraints = false, - const std::vector &starting_indices = std::vector()); + const bool use_constraints = false); /** * Computes the renumbering @@ -565,598 +430,730 @@ class DoFRenumbering * vector. */ template - static void + void compute_Cuthill_McKee (std::vector& new_dof_indices, const DH&, const bool reversed_numbering = false, - const bool use_constraints = false, - const std::vector &starting_indices = std::vector()); + const bool use_constraints = false); /** * 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 - * ommitted. + * freedom based on the BOOST + * implementation of the King + * algorithm. This often + * results in slightly larger + * (by a few percent) + * bandwidths than the + * Cuthill-McKee algorithm, + * but sparse ILUs are often + * slightly (also by a few + * percent) better + * preconditioners. * - * See the general documentation - * of this class for details on - * the different methods. - */ - template - static void - Cuthill_McKee (MGDoFHandler &dof_handler, - const unsigned int level, - const bool reversed_numbering = false, - const std::vector &starting_indices = std::vector ()); - - /** - * 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 - * component will be assigned new - * indices $n(i)$ and $n(j)$ with - * $n(i)n(j)$ if $i>j$. - * - * You can specify that the - * components are ordered in a - * different way than suggested - * 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 - * FESystem. Naming the same - * component more than once is - * possible and results in a - * blocking of several components - * into one. This is discussed in - * @ref step_22 "step-22". If you - * omit this argument, the same - * order as given by the finite - * element is used. + * As an example of the + * results of this algorithm, + * take a look at the + * comparison of various + * algorithms in the + * documentation of the + * DoFRenumbering namespace. * - * If one of the base finite - * elements from which the global - * finite element under - * consideration here, is a - * non-primitive one, i.e. its - * shape functions have more than - * one non-zero component, then - * it is not possible to - * associate these degrees of - * freedom with a single vector - * component. In this case, they - * are associated with the first - * vector component to which they - * belong. - * - * For finite elements with only - * one component, or a single - * non-primitive base element, - * this function is the identity - * operation. - */ - template - static void - component_wise (DoFHandler &dof_handler, - const std::vector &target_component - = std::vector()); - - /** - * Sort the degrees of freedom by - * component. It does the same - * thing as the above function, - * only that it does this for one - * single level of a multi-level - * discretization. The - * non-multigrid part of the - * MGDoFHandler is not touched. + * This algorithm is used in + * @ref step_22 "step-22". */ - template - static void - component_wise (MGDoFHandler& dof_handler, - unsigned int level, - const std::vector& target_component = std::vector()); - + template + void + king_ordering (DH& dof_handler, + const bool reversed_numbering = false, + const bool use_constraints = false); /** - * Sort the degrees of freedom by - * component. It does the same - * thing as the previous - * functions, but more: it - * renumbers not only every level - * of the multigrid part, but - * also the global, - * i.e. non-multigrid components. - */ - template - static void - component_wise (MGDoFHandler& dof_handler, - const std::vector& target_component = std::vector()); - - - protected: - /** - * Computes the renumbering - * vector needed by the - * component_wise() - * functions. Does not perform - * the renumbering on the - * DoFHandler dofs but returns - * the renumbering vector. - */ - template - static unsigned int - compute_component_wise (std::vector& new_dof_indices, - ITERATOR& start, - const ENDITERATOR& end, - const std::vector &target_component); - public: - /** - * Cell-wise renumbering for DG - * elements. This function takes - * the ordered set of cells in - * 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 - * Elements, i.e. all degrees of - * freedom have to be associated - * with the interior of the cell. + * Compute the renumbering + * for the King algorithm but + * do not actually renumber + * the degrees of freedom in + * the DoF handler argument. */ template - static void - cell_wise (DH &dof_handler, - const std::vector &cell_order); + void + compute_king_ordering (std::vector& new_dof_indices, + const DH&, + const bool reversed_numbering = false, + const bool use_constraints = false); /** - * @deprecated Use cell_wise() instead. + * Renumber the degrees of + * freedom based on the BOOST + * implementation of the + * minimum degree + * algorithm. Unlike the + * Cuthill-McKee algorithm, + * this algorithm does not + * attempt to minimize the + * bandwidth of a matrix but + * to minimize the amount of + * fill-in when doing an LU + * decomposition. It may + * sometimes yield better + * ILUs because of this + * property. * - * Cell-wise numbering only for DG. + * As an example of the + * results of this algorithm, + * take a look at the + * comparison of various + * algorithms in the + * documentation of the + * DoFRenumbering namespace. */ template - static void - cell_wise_dg (DH &dof_handler, - const std::vector &cell_order); + void + minimum_degree (DH& dof_handler, + const bool reversed_numbering = false, + const bool use_constraints = false); /** - * Computes the renumbering - * vector needed by the - * cell_wise_dg() function. Does - * not perform the renumbering on - * the DoFHandler dofs but - * returns the renumbering - * vector. + * Compute the renumbering + * for the minimum degree + * algorithm but do not + * actually renumber the + * degrees of freedom in the + * DoF handler argument. */ template - static void - compute_cell_wise_dg (std::vector& renumbering, - std::vector& inverse_renumbering, - const DH &dof_handler, - const std::vector &cell_order); + void + compute_minimum_degree (std::vector& new_dof_indices, + const DH&, + const bool reversed_numbering = false, + const bool use_constraints = false); + } + + /** + * 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. + * + * As an example of the results + * of this algorithm, take a look + * at the comparison of various + * algorithms in the + * documentation of the + * DoFRenumbering namespace. + */ + template + void + Cuthill_McKee (DH& dof_handler, + const bool reversed_numbering = false, + const bool use_constraints = false, + const std::vector &starting_indices = std::vector()); - /** - * Computes the renumbering - * vector needed by the - * cell_wise() function. Does - * not perform the renumbering on - * the DoFHandler dofs but - * returns the renumbering - * vector. - */ - template - static void - compute_cell_wise (std::vector& renumbering, - std::vector& inverse_renumbering, - const DH &dof_handler, - const std::vector &cell_order); + /** + * Computes the renumbering + * vector needed by the + * Cuthill_McKee() function. Does + * not perform the renumbering on + * the DoFHandler dofs but + * returns the renumbering + * vector. + */ + template + void + compute_Cuthill_McKee (std::vector& new_dof_indices, + const DH&, + const bool reversed_numbering = false, + const bool use_constraints = false, + 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. + * + * 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. + */ + template + void + Cuthill_McKee (MGDoFHandler &dof_handler, + const unsigned int level, + const bool reversed_numbering = false, + const std::vector &starting_indices = std::vector ()); + + /** + * 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 + * component will be assigned new + * indices $n(i)$ and $n(j)$ with + * $n(i)n(j)$ if $i>j$. + * + * You can specify that the + * components are ordered in a + * different way than suggested + * 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 + * FESystem. Naming the same + * component more than once is + * possible and results in a + * blocking of several components + * into one. This is discussed in + * @ref step_22 "step-22". 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 + * finite element under + * consideration here, is a + * non-primitive one, i.e. its + * shape functions have more than + * one non-zero component, then + * it is not possible to + * associate these degrees of + * freedom with a single vector + * component. In this case, they + * are associated with the first + * vector component to which they + * belong. + * + * For finite elements with only + * one component, or a single + * non-primitive base element, + * this function is the identity + * operation. + */ + template + void + component_wise (DoFHandler &dof_handler, + const std::vector &target_component + = std::vector()); + + /** + * Sort the degrees of freedom by + * component. It does the same + * thing as the above function, + * only that it does this for one + * single level of a multi-level + * discretization. The + * non-multigrid part of the + * MGDoFHandler is not touched. + */ + template + void + component_wise (MGDoFHandler& dof_handler, + unsigned int level, + const std::vector& target_component = std::vector()); + + + /** + * Sort the degrees of freedom by + * component. It does the same + * thing as the previous + * functions, but more: it + * renumbers not only every level + * of the multigrid part, but + * also the global, + * i.e. non-multigrid components. + */ + template + void + component_wise (MGDoFHandler& dof_handler, + const std::vector& target_component = std::vector()); - /** - * Cell-wise renumbering on one - * level. See the other function - * with the same name. - */ - template - static void - cell_wise (MGDoFHandler &dof_handler, - const unsigned int level, - const std::vector::cell_iterator> &cell_order); - /** - * @deprecated Use cell_wise() instead. - * - * Cell-wise renumbering on one - * level for DG elements. See the - * other function with the same - * name. - */ - template - static void - cell_wise_dg (MGDoFHandler &dof_handler, - const unsigned int level, - const std::vector::cell_iterator> &cell_order); + /** + * Computes the renumbering + * vector needed by the + * component_wise() + * functions. Does not perform + * the renumbering on the + * DoFHandler dofs but returns + * the renumbering vector. + */ + template + static unsigned int + compute_component_wise (std::vector& new_dof_indices, + ITERATOR& start, + const ENDITERATOR& end, + const std::vector &target_component); + + /** + * Cell-wise renumbering for DG + * elements. This function takes + * the ordered set of cells in + * 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 + * Elements, i.e. all degrees of + * freedom have to be associated + * with the interior of the cell. + */ + template + void + cell_wise (DH &dof_handler, + const std::vector &cell_order); + + /** + * @deprecated Use cell_wise() instead. + * + * Cell-wise numbering only for DG. + */ + template + void + cell_wise_dg (DH &dof_handler, + const std::vector &cell_order); + + /** + * Computes the renumbering + * vector needed by the + * cell_wise_dg() function. Does + * not perform the renumbering on + * the DoFHandler dofs but + * returns the renumbering + * vector. + */ + template + void + compute_cell_wise_dg (std::vector& renumbering, + std::vector& inverse_renumbering, + const DH &dof_handler, + const std::vector &cell_order); + + /** + * Computes the renumbering + * vector needed by the + * cell_wise() function. Does + * not perform the renumbering on + * the DoFHandler dofs but + * returns the renumbering + * vector. + */ + template + void + compute_cell_wise (std::vector& renumbering, + std::vector& inverse_renumbering, + const DH &dof_handler, + const std::vector &cell_order); + + /** + * Cell-wise renumbering on one + * level. See the other function + * with the same name. + */ + template + void + cell_wise (MGDoFHandler &dof_handler, + const unsigned int level, + const std::vector::cell_iterator> &cell_order); - /** - * Computes the renumbering - * vector needed by the - * cell_wise_dg() level renumbering function. Does - * not perform the renumbering on - * the MGDoFHandler dofs but - * returns the renumbering - * vector. - */ - template - static void - compute_cell_wise_dg (std::vector& renumbering, - std::vector& inverse_renumbering, - const MGDoFHandler& dof_handler, - const unsigned int level, - const std::vector::cell_iterator>& cell_order); + /** + * @deprecated Use cell_wise() instead. + * + * Cell-wise renumbering on one + * level for DG elements. See the + * other function with the same + * name. + */ + template + void + cell_wise_dg (MGDoFHandler &dof_handler, + const unsigned int level, + const std::vector::cell_iterator> &cell_order); + + /** + * Computes the renumbering + * vector needed by the + * cell_wise_dg() level renumbering function. Does + * not perform the renumbering on + * the MGDoFHandler dofs but + * returns the renumbering + * vector. + */ + template + void + compute_cell_wise_dg (std::vector& renumbering, + std::vector& inverse_renumbering, + const MGDoFHandler& dof_handler, + const unsigned int level, + const std::vector::cell_iterator>& cell_order); - /** - * Computes the renumbering - * vector needed by the - * cell_wise() level renumbering function. Does - * not perform the renumbering on - * the MGDoFHandler dofs but - * returns the renumbering - * vector. - */ - template - static void - compute_cell_wise (std::vector& renumbering, - std::vector& inverse_renumbering, - const MGDoFHandler& dof_handler, - const unsigned int level, - const std::vector::cell_iterator>& cell_order); + /** + * Computes the renumbering + * vector needed by the + * cell_wise() level renumbering function. Does + * not perform the renumbering on + * the MGDoFHandler dofs but + * returns the renumbering + * vector. + */ + template + void + compute_cell_wise (std::vector& renumbering, + std::vector& inverse_renumbering, + const MGDoFHandler& dof_handler, + const unsigned int level, + const std::vector::cell_iterator>& cell_order); - /** - * Cell-wise downstream numbering - * with respect to a constant - * flow direction. - * - * The cells are sorted such that - * the centers of higher numbers - * are further downstream with - * respect to the constant vector - * @p direction than the centers - * of lower numbers. Even if this - * yields a downstream numbering - * with respect to the flux on - * the edges for fairly general - * grids, this might not be - * guaranteed for all meshes. - * - * 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. - */ - template - static void - downstream (DH& dof_handler, - const Point& direction); + /** + * Cell-wise downstream numbering + * with respect to a constant + * flow direction. + * + * The cells are sorted such that + * the centers of higher numbers + * are further downstream with + * respect to the constant vector + * @p direction than the centers + * of lower numbers. Even if this + * yields a downstream numbering + * with respect to the flux on + * the edges for fairly general + * grids, this might not be + * guaranteed for all meshes. + * + * 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. + */ + template + void + downstream (DH& dof_handler, + const Point& direction); - /** - * @deprecated Use downstream() instead. - */ - template - static void - downstream_dg (DH& dof_handler, - const Point& direction); + /** + * @deprecated Use downstream() instead. + */ + template + void + downstream_dg (DH& dof_handler, + const Point& direction); - /** - * Cell-wise downstream numbering - * with respect to a constant - * flow direction on one - * level. See the other function - * with the same name. - */ - template - static void - downstream (MGDoFHandler &dof_handler, - const unsigned int level, - const Point &direction); - /** - * @deprecated Use downstream() - * instead. - */ - template - static void - downstream_dg (MGDoFHandler &dof_handler, - const unsigned int level, - const Point &direction); + /** + * Cell-wise downstream numbering + * with respect to a constant + * flow direction on one + * level. See the other function + * with the same name. + */ + template + void + downstream (MGDoFHandler &dof_handler, + const unsigned int level, + const Point &direction); + /** + * @deprecated Use downstream() + * instead. + */ + template + void + downstream_dg (MGDoFHandler &dof_handler, + const unsigned int level, + const Point &direction); - /** - * @deprecated The new function - * of this name computes the - * renumbering and its inverse at - * the same time. So, at least if - * you need both, you should use - * the other one. - * - * Computes the renumbering - * vector needed by the - * downstream_dg() function. Does - * not perform the renumbering on - * the DoFHandler dofs but - * returns the renumbering - * vector. - */ - template - static void - compute_downstream_dg (std::vector& new_dof_indices, - const DH& dof_handler, - const Point& direction); + /** + * @deprecated The new function + * of this name computes the + * renumbering and its inverse at + * the same time. So, at least if + * you need both, you should use + * the other one. + * + * Computes the renumbering + * vector needed by the + * downstream_dg() function. Does + * not perform the renumbering on + * the DoFHandler dofs but + * returns the renumbering + * vector. + */ + template + void + compute_downstream_dg (std::vector& new_dof_indices, + const DH& dof_handler, + const Point& direction); - /** - * Computes the renumbering - * vector needed by the - * downstream_dg() function. Does - * not perform the renumbering on - * the DoFHandler dofs but - * returns the renumbering - * vector. - */ - template - static void - compute_downstream (std::vector& new_dof_indices, - std::vector& reverse, - const DH& dof_handler, - const Point& direction); + /** + * Computes the renumbering + * vector needed by the + * downstream_dg() function. Does + * not perform the renumbering on + * the DoFHandler dofs but + * returns the renumbering + * vector. + */ + template + void + compute_downstream (std::vector& new_dof_indices, + std::vector& reverse, + const DH& dof_handler, + const Point& direction); - /** - * @deprecated Use - * compute_downstream() instead - */ - template - static void - compute_downstream_dg (std::vector& new_dof_indices, - std::vector& reverse, - const DH& dof_handler, - const Point& direction); + /** + * @deprecated Use + * compute_downstream() instead + */ + template + void + compute_downstream_dg (std::vector& new_dof_indices, + std::vector& reverse, + const DH& dof_handler, + const Point& direction); - /** - * Computes the renumbering - * vector needed by the - * downstream_dg() function. Does - * not perform the renumbering on - * the MGDoFHandler dofs but - * returns the renumbering - * vector. - */ - template - static void - compute_downstream (std::vector& new_dof_indices, - std::vector& reverse, - const MGDoFHandler& dof_handler, - const unsigned int level, - const Point& direction); + /** + * Computes the renumbering + * vector needed by the + * downstream_dg() function. Does + * not perform the renumbering on + * the MGDoFHandler dofs but + * returns the renumbering + * vector. + */ + template + void + compute_downstream (std::vector& new_dof_indices, + std::vector& reverse, + const MGDoFHandler& dof_handler, + const unsigned int level, + const Point& direction); - /** - * @deprecated Use - * compute_downstream() instead - */ - template - static void - compute_downstream_dg (std::vector& new_dof_indices, - std::vector& reverse, - const MGDoFHandler& dof_handler, - const unsigned int level, - const Point& direction); + /** + * @deprecated Use + * compute_downstream() instead + */ + template + void + compute_downstream_dg (std::vector& new_dof_indices, + std::vector& reverse, + const MGDoFHandler& dof_handler, + const unsigned int level, + const Point& direction); - /** - * Cell-wise clockwise numbering. - * - * This function produces a - * (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. - */ - template - static void - clockwise_dg (DH& dof_handler, - const Point& center, - const bool counter = false); + /** + * Cell-wise clockwise numbering. + * + * This function produces a + * (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. + */ + template + void + clockwise_dg (DH& dof_handler, + const Point& center, + const bool counter = false); - /** - * Cell-wise clockwise numbering - * on one level. See the other - * function with the same name. - */ - template - static void - clockwise_dg (MGDoFHandler &dof_handler, - const unsigned int level, - const Point ¢er, - const bool counter = false); + /** + * Cell-wise clockwise numbering + * on one level. See the other + * function with the same name. + */ + template + void + clockwise_dg (MGDoFHandler &dof_handler, + const unsigned int level, + const Point ¢er, + const bool counter = false); - /** - * Computes the renumbering - * vector needed by the - * clockwise_dg() functions. Does - * not perform the renumbering on - * the DoFHandler dofs but - * returns the renumbering - * vector. - */ - template - static void - compute_clockwise_dg (std::vector& new_dof_indices, - const DH& dof_handler, - const Point& center, - const bool counter); + /** + * Computes the renumbering + * vector needed by the + * clockwise_dg() functions. Does + * not perform the renumbering on + * the DoFHandler dofs but + * returns the renumbering + * vector. + */ + template + void + compute_clockwise_dg (std::vector& new_dof_indices, + const DH& dof_handler, + const Point& center, + const bool counter); - /** - * Sort those degrees of freedom - * which are tagged with @p true - * in the @p selected_dofs array - * to the back of the DoF - * numbers. The sorting is - * stable, i.e. the relative - * order within the tagged - * degrees of freedom is - * preserved, as is the relative - * order within the untagged - * ones. - */ - template - static void - sort_selected_dofs_back (DH& dof_handler, - const std::vector& selected_dofs); + /** + * Sort those degrees of freedom + * which are tagged with @p true + * in the @p selected_dofs array + * to the back of the DoF + * numbers. The sorting is + * stable, i.e. the relative + * order within the tagged + * degrees of freedom is + * preserved, as is the relative + * order within the untagged + * ones. + */ + template + void + sort_selected_dofs_back (DH& dof_handler, + const std::vector& selected_dofs); - /** - * Computes the renumbering - * vector needed by the - * sort_selected_dofs_back() - * function. Does not perform the - * renumbering on the DoFHandler - * dofs but returns the - * renumbering vector. - */ - template - static void - compute_sort_selected_dofs_back (std::vector& new_dof_indices, - const DH& dof_handler, - const std::vector& selected_dofs); + /** + * Computes the renumbering + * vector needed by the + * sort_selected_dofs_back() + * function. Does not perform the + * renumbering on the DoFHandler + * dofs but returns the + * renumbering vector. + */ + template + void + compute_sort_selected_dofs_back (std::vector& new_dof_indices, + const DH& dof_handler, + const std::vector& selected_dofs); - /** - * Renumber the degrees of - * freedom in a random way. - */ - template - static void - random (DH& dof_handler); + /** + * Renumber the degrees of + * freedom in a random way. + */ + template + void + random (DH& dof_handler); - /** - * Computes the renumbering - * vector needed by the random() - * function. Does not perform the - * renumbering on the DoFHandler - * dofs but returns the - * renumbering vector. - */ - template - static void - compute_random (std::vector &new_dof_indices, - const DH& dof_handler); + /** + * Computes the renumbering + * vector needed by the random() + * function. Does not perform the + * renumbering on the DoFHandler + * dofs but returns the + * renumbering vector. + */ + template + void + compute_random (std::vector &new_dof_indices, + const DH& 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 - * 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 - * DoFTools::get_subdomain_association - * function. - * - * The algorithm is stable, i.e. if - * two dofs i,j have i and belong - * to the same subdomain, then they - * will be in this order also after - * reordering. - */ - template - static void - subdomain_wise (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 + * 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 + * DoFTools::get_subdomain_association + * function. + * + * The algorithm is stable, i.e. if + * two dofs i,j have i and belong + * to the same subdomain, then they + * will be in this order also after + * reordering. + */ + template + void + subdomain_wise (DoFHandler &dof_handler); - /** - * 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. - */ - template - static void - compute_subdomain_wise (std::vector &new_dof_indices, - const DoFHandler &dof_handler); + /** + * 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. + */ + template + void + compute_subdomain_wise (std::vector &new_dof_indices, + const DoFHandler &dof_handler); - /** - * Exception - * - * @ingroup Exceptions - */ - DeclException0 (ExcRenumberingIncomplete); - /** - * Exception - * - * @ingroup Exceptions - */ - DeclException0 (ExcInvalidComponentOrder); - /** - * The function is only - * implemented for Discontinuous - * Galerkin Finite elements. - * - * @ingroup Exceptions - */ - DeclException0 (ExcNotDGFEM); -}; + /** + * Exception + * + * @ingroup Exceptions + */ + DeclException0 (ExcRenumberingIncomplete); + /** + * Exception + * + * @ingroup Exceptions + */ + DeclException0 (ExcInvalidComponentOrder); + /** + * The function is only + * implemented for Discontinuous + * Galerkin Finite elements. + * + * @ingroup Exceptions + */ + DeclException0 (ExcNotDGFEM); +} DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/deal.II/source/dofs/dof_renumbering.cc b/deal.II/deal.II/source/dofs/dof_renumbering.cc index b0e1eebd95..b38a596eab 100644 --- a/deal.II/deal.II/source/dofs/dof_renumbering.cc +++ b/deal.II/deal.II/source/dofs/dof_renumbering.cc @@ -63,55 +63,83 @@ extern "C" long int lrand48 (void); DEAL_II_NAMESPACE_OPEN +namespace internal +{ + // the following function is sorta + // stupid: gcc complains that the + // variables _3, _4, ... declared + // in + // contrib/boost/include/boost/bind/placeholders.hpp + // are unused. _1 and _2 are + // apparently used, so this isn't a + // problem. to suppress the warning + // in this case, write a dummy + // function that simply references + // these variables in an + // unobtrusive way... + void avoid_warning_about_unused_variables () + { + std::printf ("%p %p %p %p %p %p %p", + (void*)&_3, (void*)&_4, + (void*)&_5,(void*) &_6, + (void*)&_7,(void*) &_8, + (void*)&_9); + } +} + +namespace DoFRenumbering +{ + namespace internal + { // The following two classes are defined to be used in the compute_* // functions. Using them allows to use the same function to compute // level numbering for multigrid as well as numbering for the global // system matrix. -template -class WrapDoFIterator : private T -{ - public: - WrapDoFIterator (const T& t) : T(t) {} - - void get_dof_indices (std::vector& v) const - { - (*this)->get_dof_indices(v); - } - - template - bool operator != (const T2& i) const - { - return (! (T::operator==(i))); - } - // Allow access to these operators of T - T::operator->; - T::operator++; - T::operator==; -}; - - + template + class WrapDoFIterator : private T + { + public: + WrapDoFIterator (const T& t) : T(t) {} + + void get_dof_indices (std::vector& v) const + { + (*this)->get_dof_indices(v); + } + + template + bool operator != (const T2& i) const + { + return (! (T::operator==(i))); + } + // Allow access to these operators of T + T::operator->; + T::operator++; + T::operator==; + }; -template -class WrapMGDoFIterator : private T -{ - public: - WrapMGDoFIterator (const T& t) : T(t) {} - void get_dof_indices (std::vector& v) const - { - (*this)->get_mg_dof_indices(v); - } - bool operator != (const WrapMGDoFIterator& i) const - { - return (! (T::operator==(i))); - } - // Allow access to these operators of T - T::operator->; - T::operator++; - T::operator==; -}; + template + class WrapMGDoFIterator : private T + { + public: + WrapMGDoFIterator (const T& t) : T(t) {} + + void get_dof_indices (std::vector& v) const + { + (*this)->get_mg_dof_indices(v); + } + + bool operator != (const WrapMGDoFIterator& i) const + { + return (! (T::operator==(i))); + } + // Allow access to these operators of T + T::operator->; + T::operator++; + T::operator==; + }; @@ -120,2016 +148,2032 @@ class WrapMGDoFIterator : private T * of the second cell is downstream of the center of the first one with * respect to the direction given to the constructor. */ -template -struct CompareDownstream -{ - /** - * Constructor. - */ - CompareDownstream (const Point &dir) - : - dir(dir) - {} - /** - * Return true if c1 < c2. - */ - template - bool operator () (const DHCellIterator &c1, - const DHCellIterator &c2) const - { - const Point diff = c2->center() - c1->center(); - return (diff*dir > 0); - } - - private: - /** - * Flow direction. - */ - const Point dir; -}; - - -namespace types -{ - using namespace boost; - using namespace std; + template + struct CompareDownstream + { + /** + * Constructor. + */ + CompareDownstream (const Point &dir) + : + dir(dir) + {} + /** + * Return true if c1 < c2. + */ + template + bool operator () (const DHCellIterator &c1, + const DHCellIterator &c2) const + { + const Point diff = c2->center() - c1->center(); + return (diff*dir > 0); + } + + private: + /** + * Flow direction. + */ + const Point dir; + }; + } + + namespace boost + { + namespace types + { + using namespace ::boost; + using namespace std; - typedef vecS x; + typedef vecS x; - typedef adjacency_list > > Graph; - typedef graph_traits::vertex_descriptor Vertex; - typedef graph_traits::vertices_size_type size_type; + typedef adjacency_list > > Graph; + typedef graph_traits::vertex_descriptor Vertex; + typedef graph_traits::vertices_size_type size_type; - typedef std::pair Pair; -} + typedef std::pair Pair; + } -namespace -{ - template - void create_graph (const DH &dof_handler, - const bool use_constraints, - types::Graph &graph, - types::property_map::type &graph_degree) - { - Assert (use_constraints == false, ExcNotImplemented()); - - std::vector dofs_on_this_cell; - dofs_on_this_cell.reserve (DoFTools::max_dofs_per_cell(dof_handler)); - typename DH::active_cell_iterator cell = dof_handler.begin_active(), - endc = dof_handler.end(); - for (; cell!=endc; ++cell) + namespace internal + { + template + void create_graph (const DH &dof_handler, + const bool use_constraints, + types::Graph &graph, + types::property_map::type &graph_degree) { - const unsigned int dofs_per_cell = cell->get_fe().dofs_per_cell; - dofs_on_this_cell.resize (dofs_per_cell); - cell->get_dof_indices (dofs_on_this_cell); + Assert (use_constraints == false, ExcNotImplemented()); + + std::vector dofs_on_this_cell; + dofs_on_this_cell.reserve (DoFTools::max_dofs_per_cell(dof_handler)); + typename DH::active_cell_iterator cell = dof_handler.begin_active(), + endc = dof_handler.end(); + for (; cell!=endc; ++cell) + { + const unsigned int dofs_per_cell = cell->get_fe().dofs_per_cell; + dofs_on_this_cell.resize (dofs_per_cell); + cell->get_dof_indices (dofs_on_this_cell); - for (unsigned int i=0; i::vertex_iterator ui, ui_end; + types::graph_traits::vertex_iterator ui, ui_end; - graph_degree = get(boost::vertex_degree, graph); - for (boost::tie(ui, ui_end) = vertices(graph); ui != ui_end; ++ui) - graph_degree[*ui] = degree(*ui, graph); + graph_degree = get(::boost::vertex_degree, graph); + for (::boost::tie(ui, ui_end) = vertices(graph); ui != ui_end; ++ui) + graph_degree[*ui] = degree(*ui, graph); - } -} + } + } -template -void -DoFRenumbering::boost:: -Cuthill_McKee (DH& dof_handler, - const bool reversed_numbering, - const bool use_constraints) -{ - std::vector renumbering(dof_handler.n_dofs(), - DH::invalid_dof_index); - compute_Cuthill_McKee(renumbering, dof_handler, reversed_numbering, - use_constraints); - - // actually perform renumbering; - // this is dimension specific and - // thus needs an own function - dof_handler.renumber_dofs (renumbering); -} + template + void + Cuthill_McKee (DH& dof_handler, + const bool reversed_numbering, + const bool use_constraints) + { + std::vector renumbering(dof_handler.n_dofs(), + DH::invalid_dof_index); + compute_Cuthill_McKee(renumbering, dof_handler, reversed_numbering, + use_constraints); + + // actually perform renumbering; + // this is dimension specific and + // thus needs an own function + dof_handler.renumber_dofs (renumbering); + } -template -void -DoFRenumbering::boost:: -compute_Cuthill_McKee (std::vector& new_dof_indices, - const DH &dof_handler, - const bool reversed_numbering, - const bool use_constraints) -{ - types::Graph - graph(dof_handler.n_dofs()); - types::property_map::type - graph_degree; + template + void + compute_Cuthill_McKee (std::vector& new_dof_indices, + const DH &dof_handler, + const bool reversed_numbering, + const bool use_constraints) + { + types::Graph + graph(dof_handler.n_dofs()); + types::property_map::type + graph_degree; - create_graph (dof_handler, use_constraints, graph, graph_degree); + internal::create_graph (dof_handler, use_constraints, graph, graph_degree); - types::property_map::type - index_map = get(::boost::vertex_index, graph); + types::property_map::type + index_map = get(::boost::vertex_index, graph); - std::vector inv_perm(num_vertices(graph)); + std::vector inv_perm(num_vertices(graph)); - if (reversed_numbering == false) - ::boost::cuthill_mckee_ordering(graph, inv_perm.rbegin(), - get(::boost::vertex_color, graph), - make_degree_map(graph)); - else - ::boost::cuthill_mckee_ordering(graph, inv_perm.begin(), - get(::boost::vertex_color, graph), - make_degree_map(graph)); + if (reversed_numbering == false) + ::boost::cuthill_mckee_ordering(graph, inv_perm.rbegin(), + get(::boost::vertex_color, graph), + make_degree_map(graph)); + else + ::boost::cuthill_mckee_ordering(graph, inv_perm.begin(), + get(::boost::vertex_color, graph), + make_degree_map(graph)); - for (types::size_type c = 0; c != inv_perm.size(); ++c) - new_dof_indices[index_map[inv_perm[c]]] = c; + for (types::size_type c = 0; c != inv_perm.size(); ++c) + new_dof_indices[index_map[inv_perm[c]]] = c; - Assert (std::find (new_dof_indices.begin(), new_dof_indices.end(), - DH::invalid_dof_index) == new_dof_indices.end(), - ExcInternalError()); -} + Assert (std::find (new_dof_indices.begin(), new_dof_indices.end(), + DH::invalid_dof_index) == new_dof_indices.end(), + ExcInternalError()); + } -template -void -DoFRenumbering::boost:: -king_ordering (DH& dof_handler, - const bool reversed_numbering, - const bool use_constraints) -{ - std::vector renumbering(dof_handler.n_dofs(), - DH::invalid_dof_index); - compute_king_ordering(renumbering, dof_handler, reversed_numbering, - use_constraints); - - // actually perform renumbering; - // this is dimension specific and - // thus needs an own function - dof_handler.renumber_dofs (renumbering); -} + template + void + king_ordering (DH& dof_handler, + const bool reversed_numbering, + const bool use_constraints) + { + std::vector renumbering(dof_handler.n_dofs(), + DH::invalid_dof_index); + compute_king_ordering(renumbering, dof_handler, reversed_numbering, + use_constraints); + + // actually perform renumbering; + // this is dimension specific and + // thus needs an own function + dof_handler.renumber_dofs (renumbering); + } -template -void -DoFRenumbering::boost:: -compute_king_ordering (std::vector& new_dof_indices, - const DH &dof_handler, - const bool reversed_numbering, - const bool use_constraints) -{ - types::Graph - graph(dof_handler.n_dofs()); - types::property_map::type - graph_degree; + template + void + compute_king_ordering (std::vector& new_dof_indices, + const DH &dof_handler, + const bool reversed_numbering, + const bool use_constraints) + { + types::Graph + graph(dof_handler.n_dofs()); + types::property_map::type + graph_degree; - create_graph (dof_handler, use_constraints, graph, graph_degree); + internal::create_graph (dof_handler, use_constraints, graph, graph_degree); - types::property_map::type - index_map = get(::boost::vertex_index, graph); + types::property_map::type + index_map = get(::boost::vertex_index, graph); - std::vector inv_perm(num_vertices(graph)); + std::vector inv_perm(num_vertices(graph)); - if (reversed_numbering == false) - ::boost::king_ordering(graph, inv_perm.rbegin()); - else - ::boost::king_ordering(graph, inv_perm.begin()); + if (reversed_numbering == false) + ::boost::king_ordering(graph, inv_perm.rbegin()); + else + ::boost::king_ordering(graph, inv_perm.begin()); - for (types::size_type c = 0; c != inv_perm.size(); ++c) - new_dof_indices[index_map[inv_perm[c]]] = c; + for (types::size_type c = 0; c != inv_perm.size(); ++c) + new_dof_indices[index_map[inv_perm[c]]] = c; - Assert (std::find (new_dof_indices.begin(), new_dof_indices.end(), - DH::invalid_dof_index) == new_dof_indices.end(), - ExcInternalError()); -} + Assert (std::find (new_dof_indices.begin(), new_dof_indices.end(), + DH::invalid_dof_index) == new_dof_indices.end(), + ExcInternalError()); + } -template -void -DoFRenumbering::boost:: -minimum_degree (DH& dof_handler, - const bool reversed_numbering, - const bool use_constraints) -{ - std::vector renumbering(dof_handler.n_dofs(), - DH::invalid_dof_index); - compute_minimum_degree(renumbering, dof_handler, reversed_numbering, - use_constraints); - - // actually perform renumbering; - // this is dimension specific and - // thus needs an own function - dof_handler.renumber_dofs (renumbering); -} + template + void + minimum_degree (DH& dof_handler, + const bool reversed_numbering, + const bool use_constraints) + { + std::vector renumbering(dof_handler.n_dofs(), + DH::invalid_dof_index); + compute_minimum_degree(renumbering, dof_handler, reversed_numbering, + use_constraints); + + // actually perform renumbering; + // this is dimension specific and + // thus needs an own function + dof_handler.renumber_dofs (renumbering); + } -template -void -DoFRenumbering::boost:: -compute_minimum_degree (std::vector& new_dof_indices, - const DH &dof_handler, - const bool reversed_numbering, - const bool use_constraints) -{ - Assert (use_constraints == false, ExcNotImplemented()); + template + void + compute_minimum_degree (std::vector& new_dof_indices, + const DH &dof_handler, + const bool reversed_numbering, + const bool use_constraints) + { + Assert (use_constraints == false, ExcNotImplemented()); - // the following code is pretty - // much a verbatim copy of the - // sample code for the - // minimum_degree_ordering manual - // page from the BOOST Graph - // Library - using namespace boost; + // the following code is pretty + // much a verbatim copy of the + // sample code for the + // minimum_degree_ordering manual + // page from the BOOST Graph + // Library + using namespace ::boost; - int delta = 0; + int delta = 0; - typedef double Type; + typedef double Type; - // must be BGL directed graph now - typedef adjacency_list Graph; - typedef graph_traits::vertex_descriptor Vertex; + // must be BGL directed graph now + typedef adjacency_list Graph; + typedef graph_traits::vertex_descriptor Vertex; - int n = dof_handler.n_dofs(); + int n = dof_handler.n_dofs(); - Graph G(n); + Graph G(n); - std::vector dofs_on_this_cell; + std::vector dofs_on_this_cell; - typename DH::active_cell_iterator cell = dof_handler.begin_active(), - endc = dof_handler.end(); + typename DH::active_cell_iterator cell = dof_handler.begin_active(), + endc = dof_handler.end(); - for (; cell!=endc; ++cell) - { + for (; cell!=endc; ++cell) + { - const unsigned int dofs_per_cell = cell->get_fe().dofs_per_cell; + const unsigned int dofs_per_cell = cell->get_fe().dofs_per_cell; - dofs_on_this_cell.resize (dofs_per_cell); + dofs_on_this_cell.resize (dofs_per_cell); - cell->get_dof_indices (dofs_on_this_cell); - for (unsigned int i=0; i dofs_on_this_cell[j]) - { - add_edge (dofs_on_this_cell[i], dofs_on_this_cell[j], G); - add_edge (dofs_on_this_cell[j], dofs_on_this_cell[i], G); - } - } + cell->get_dof_indices (dofs_on_this_cell); + for (unsigned int i=0; i dofs_on_this_cell[j]) + { + add_edge (dofs_on_this_cell[i], dofs_on_this_cell[j], G); + add_edge (dofs_on_this_cell[j], dofs_on_this_cell[i], G); + } + } - typedef std::vector Vector; + typedef std::vector Vector; - Vector inverse_perm(n, 0); + Vector inverse_perm(n, 0); - Vector perm(n, 0); + Vector perm(n, 0); - Vector supernode_sizes(n, 1); - // init has to be 1 + Vector supernode_sizes(n, 1); + // init has to be 1 - ::boost::property_map::type - id = get(vertex_index, G); + ::boost::property_map::type + id = get(vertex_index, G); - Vector degree(n, 0); + Vector degree(n, 0); - minimum_degree_ordering - (G, - make_iterator_property_map(°ree[0], id, degree[0]), - &inverse_perm[0], - &perm[0], - make_iterator_property_map(&supernode_sizes[0], id, supernode_sizes[0]), - delta, id); + minimum_degree_ordering + (G, + make_iterator_property_map(°ree[0], id, degree[0]), + &inverse_perm[0], + &perm[0], + make_iterator_property_map(&supernode_sizes[0], id, supernode_sizes[0]), + delta, id); - for (int i=0; i -void -DoFRenumbering:: -Cuthill_McKee (DH& dof_handler, - const bool reversed_numbering, - const bool use_constraints, - const std::vector &starting_indices) -{ - std::vector renumbering(dof_handler.n_dofs(), - DH::invalid_dof_index); - compute_Cuthill_McKee(renumbering, dof_handler, reversed_numbering, - use_constraints, starting_indices); - - // actually perform renumbering; - // this is dimension specific and - // thus needs an own function - dof_handler.renumber_dofs (renumbering); -} + template + void + Cuthill_McKee (DH& dof_handler, + const bool reversed_numbering, + const bool use_constraints, + const std::vector &starting_indices) + { + std::vector renumbering(dof_handler.n_dofs(), + DH::invalid_dof_index); + compute_Cuthill_McKee(renumbering, dof_handler, reversed_numbering, + use_constraints, starting_indices); + + // actually perform renumbering; + // this is dimension specific and + // thus needs an own function + dof_handler.renumber_dofs (renumbering); + } -template -void -DoFRenumbering:: -compute_Cuthill_McKee (std::vector& new_indices, - const DH& dof_handler, - const bool reversed_numbering, - const bool use_constraints, - const std::vector& starting_indices) -{ - // make the connection graph. in - // more than 2d use an intermediate - // compressed sparsity pattern - // since the we don't have very - // good estimates for - // max_couplings_between_dofs() in - // 3d and this then leads to - // excessive memory consumption - // - // note that if constraints are not - // requested, then the - // 'constraints' object will be - // empty, and calling condense with - // it is a no-op - ConstraintMatrix constraints; - if (use_constraints) - DoFTools::make_hanging_node_constraints (dof_handler, constraints); - constraints.close (); + template + void + compute_Cuthill_McKee (std::vector& new_indices, + const DH& dof_handler, + const bool reversed_numbering, + const bool use_constraints, + const std::vector& starting_indices) + { + // make the connection graph. in + // more than 2d use an intermediate + // compressed sparsity pattern + // since the we don't have very + // good estimates for + // max_couplings_between_dofs() in + // 3d and this then leads to + // excessive memory consumption + // + // note that if constraints are not + // requested, then the + // 'constraints' object will be + // empty, and calling condense with + // it is a no-op + ConstraintMatrix constraints; + if (use_constraints) + DoFTools::make_hanging_node_constraints (dof_handler, constraints); + constraints.close (); - SparsityPattern sparsity; - if (DH::dimension <= 2) - { - sparsity.reinit (dof_handler.n_dofs(), - dof_handler.n_dofs(), - dof_handler.max_couplings_between_dofs()); - DoFTools::make_sparsity_pattern (dof_handler, sparsity); - constraints.condense (sparsity); - } - else - { - CompressedSparsityPattern csp (dof_handler.n_dofs(), - dof_handler.n_dofs()); - DoFTools::make_sparsity_pattern (dof_handler, csp); - constraints.condense (csp); - sparsity.copy_from (csp); - } + SparsityPattern sparsity; + if (DH::dimension <= 2) + { + sparsity.reinit (dof_handler.n_dofs(), + dof_handler.n_dofs(), + dof_handler.max_couplings_between_dofs()); + DoFTools::make_sparsity_pattern (dof_handler, sparsity); + constraints.condense (sparsity); + } + else + { + CompressedSparsityPattern csp (dof_handler.n_dofs(), + dof_handler.n_dofs()); + DoFTools::make_sparsity_pattern (dof_handler, csp); + constraints.condense (csp); + sparsity.copy_from (csp); + } - // constraints are not needed anymore - constraints.clear (); + // constraints are not needed anymore + constraints.clear (); - const unsigned int n_dofs = sparsity.n_rows(); - // store the new dof numbers; - // invalid_dof_index means that no new - // number was chosen yet - Assert(new_indices.size() == n_dofs, - ExcDimensionMismatch(new_indices.size(), n_dofs)); + const unsigned int n_dofs = sparsity.n_rows(); + // store the new dof numbers; + // invalid_dof_index means that no new + // number was chosen yet + Assert(new_indices.size() == n_dofs, + ExcDimensionMismatch(new_indices.size(), n_dofs)); - // store the indices of the dofs renumbered - // in the last round. Default to starting - // points - std::vector last_round_dofs (starting_indices); + // store the indices of the dofs renumbered + // in the last round. Default to starting + // points + std::vector last_round_dofs (starting_indices); - // delete disallowed elements - for (unsigned int i=0; i=n_dofs)) - last_round_dofs[i] = DH::invalid_dof_index; + // delete disallowed elements + for (unsigned int i=0; i=n_dofs)) + last_round_dofs[i] = DH::invalid_dof_index; - std::remove_if (last_round_dofs.begin(), last_round_dofs.end(), - std::bind2nd(std::equal_to(), - DH::invalid_dof_index)); + std::remove_if (last_round_dofs.begin(), last_round_dofs.end(), + std::bind2nd(std::equal_to(), + DH::invalid_dof_index)); - // now if no valid points remain: - // find dof with lowest coordination - // number + // now if no valid points remain: + // find dof with lowest coordination + // number - if (last_round_dofs.size() == 0) - { - unsigned int starting_point = DH::invalid_dof_index; - unsigned int min_coordination = n_dofs; - for (unsigned int row=0; row next_round_dofs; - - // find all neighbors of the - // dofs numbered in the last - // round - for (unsigned int i=0; i::iterator end_sorted; - end_sorted = std::unique (next_round_dofs.begin(), next_round_dofs.end()); - next_round_dofs.erase (end_sorted, next_round_dofs.end()); - - // eliminate dofs which are - // already numbered - for (int s=next_round_dofs.size()-1; s>=0; --s) - if (new_indices[next_round_dofs[s]] != DH::invalid_dof_index) - next_round_dofs.erase (next_round_dofs.begin() + s); - - // check whether there are any new - // dofs - all_dofs_renumbered = (next_round_dofs.size() == 0); - if (all_dofs_renumbered) - // end loop if possible - continue; - - - // store for each coordination - // number the dofs with these - // coordination number - std::multimap dofs_by_coordination; - - // find coordination number for - // each of these dofs - for (std::vector::iterator s=next_round_dofs.begin(); - s!=next_round_dofs.end(); ++s) - { - unsigned int coordination = 0; - for (unsigned int j=sparsity.get_rowstart_indices()[*s]; - j next_round_dofs; + + // find all neighbors of the + // dofs numbered in the last + // round + for (unsigned int i=0; i new_entry (coordination, *s); - dofs_by_coordination.insert (new_entry); - }; + next_round_dofs.push_back (sparsity.get_column_numbers()[j]); - // assign new DoF numbers to - // the elements of the present - // front: - std::multimap::iterator i; - for (i = dofs_by_coordination.begin(); i!=dofs_by_coordination.end(); ++i) - new_indices[i->second] = next_free_number++; - - // after that: copy this round's - // dofs for the next round - last_round_dofs = next_round_dofs; - }; + // sort dof numbers + std::sort (next_round_dofs.begin(), next_round_dofs.end()); + + // delete multiple entries + std::vector::iterator end_sorted; + end_sorted = std::unique (next_round_dofs.begin(), next_round_dofs.end()); + next_round_dofs.erase (end_sorted, next_round_dofs.end()); + + // eliminate dofs which are + // already numbered + for (int s=next_round_dofs.size()-1; s>=0; --s) + if (new_indices[next_round_dofs[s]] != DH::invalid_dof_index) + next_round_dofs.erase (next_round_dofs.begin() + s); + + // check whether there are any new + // dofs + all_dofs_renumbered = (next_round_dofs.size() == 0); + if (all_dofs_renumbered) + // end loop if possible + continue; + + + // store for each coordination + // number the dofs with these + // coordination number + std::multimap dofs_by_coordination; + + // find coordination number for + // each of these dofs + for (std::vector::iterator s=next_round_dofs.begin(); + s!=next_round_dofs.end(); ++s) + { + unsigned int coordination = 0; + for (unsigned int j=sparsity.get_rowstart_indices()[*s]; + j new_entry (coordination, *s); + dofs_by_coordination.insert (new_entry); + }; + + // assign new DoF numbers to + // the elements of the present + // front: + std::multimap::iterator i; + for (i = dofs_by_coordination.begin(); i!=dofs_by_coordination.end(); ++i) + new_indices[i->second] = next_free_number++; + + // after that: copy this round's + // dofs for the next round + last_round_dofs = next_round_dofs; + }; #ifdef DEBUG - // test for all indices - // numbered. this mostly tests - // whether the - // front-marching-algorithm (which - // Cuthill-McKee actually is) has - // reached all points. it should - // usually do so, but might not for - // two reasons: - // - // - The algorithm above has a bug, or - // - The domain is not connected and - // consists of separate parts. - // - // In any case, if not all DoFs - // have been reached, renumbering - // will not be possible - if (std::find (new_indices.begin(), new_indices.end(), DH::invalid_dof_index) - != - new_indices.end()) - Assert (false, ExcRenumberingIncomplete()); - Assert (next_free_number == n_dofs, - ExcRenumberingIncomplete()); + // test for all indices + // numbered. this mostly tests + // whether the + // front-marching-algorithm (which + // Cuthill-McKee actually is) has + // reached all points. it should + // usually do so, but might not for + // two reasons: + // + // - The algorithm above has a bug, or + // - The domain is not connected and + // consists of separate parts. + // + // In any case, if not all DoFs + // have been reached, renumbering + // will not be possible + if (std::find (new_indices.begin(), new_indices.end(), DH::invalid_dof_index) + != + new_indices.end()) + Assert (false, ExcRenumberingIncomplete()); + Assert (next_free_number == n_dofs, + ExcRenumberingIncomplete()); #endif - if (reversed_numbering) - for (std::vector::iterator i=new_indices.begin(); - i!=new_indices.end(); ++i) - *i = n_dofs-*i-1; -} + if (reversed_numbering) + for (std::vector::iterator i=new_indices.begin(); + i!=new_indices.end(); ++i) + *i = n_dofs-*i-1; + } -template -void DoFRenumbering::Cuthill_McKee ( - MGDoFHandler &dof_handler, - const unsigned int level, - const bool reversed_numbering, - const std::vector &starting_indices) -{ + template + void Cuthill_McKee (MGDoFHandler &dof_handler, + const unsigned int level, + const bool reversed_numbering, + const std::vector &starting_indices) + { //TODO: we should be doing the same here as in the other compute_CMK function to preserve some memory - // make the connection graph - SparsityPattern sparsity (dof_handler.n_dofs(level), - dof_handler.max_couplings_between_dofs()); - MGTools::make_sparsity_pattern (dof_handler, sparsity, level); + // make the connection graph + SparsityPattern sparsity (dof_handler.n_dofs(level), + dof_handler.max_couplings_between_dofs()); + MGTools::make_sparsity_pattern (dof_handler, sparsity, level); - const unsigned int n_dofs = sparsity.n_rows(); - // store the new dof numbers; invalid_dof_index means - // that no new number was chosen yet - std::vector new_indices(n_dofs, DoFHandler::invalid_dof_index); - - // store the indices of the dofs renumbered - // in the last round. Default to starting - // points - std::vector last_round_dofs (starting_indices); - - // delete disallowed elements - for (unsigned int i=0; i::invalid_dof_index) || - (last_round_dofs[i]>=n_dofs)) - last_round_dofs[i] = DoFHandler::invalid_dof_index; - - std::remove_if (last_round_dofs.begin(), last_round_dofs.end(), - std::bind2nd(std::equal_to(), - DoFHandler::invalid_dof_index)); - - // now if no valid points remain: - // find dof with lowest coordination - // number - - if (last_round_dofs.size() == 0) - { - unsigned int starting_point = DoFHandler::invalid_dof_index; - unsigned int min_coordination = n_dofs; - for (unsigned int row=0; row new_indices(n_dofs, DoFHandler::invalid_dof_index); + + // store the indices of the dofs renumbered + // in the last round. Default to starting + // points + std::vector last_round_dofs (starting_indices); + + // delete disallowed elements + for (unsigned int i=0; i::invalid_dof_index) || + (last_round_dofs[i]>=n_dofs)) + last_round_dofs[i] = DoFHandler::invalid_dof_index; + + std::remove_if (last_round_dofs.begin(), last_round_dofs.end(), + std::bind2nd(std::equal_to(), + DoFHandler::invalid_dof_index)); + + // now if no valid points remain: + // find dof with lowest coordination + // number + + if (last_round_dofs.size() == 0) + { + unsigned int starting_point = DoFHandler::invalid_dof_index; + unsigned int min_coordination = n_dofs; + for (unsigned int row=0; row next_round_dofs; - - // find all neighbors of the - // dofs numbered in the last - // round - for (unsigned int i=0; i::iterator end_sorted; - end_sorted = std::unique (next_round_dofs.begin(), next_round_dofs.end()); - next_round_dofs.erase (end_sorted, next_round_dofs.end()); - - // eliminate dofs which are - // already numbered - for (int s=next_round_dofs.size()-1; s>=0; --s) - if (new_indices[next_round_dofs[s]] != DoFHandler::invalid_dof_index) - next_round_dofs.erase (next_round_dofs.begin() + s); - - // check whether there are any new - // dofs - all_dofs_renumbered = (next_round_dofs.size() == 0); - if (all_dofs_renumbered) - // end loop if possible - continue; - - - // store for each coordination - // number the dofs with these - // coordination number - std::multimap dofs_by_coordination; - - // find coordination number for - // each of these dofs - for (std::vector::iterator s=next_round_dofs.begin(); - s!=next_round_dofs.end(); ++s) - { - unsigned int coordination = 0; - for (unsigned int j=sparsity.get_rowstart_indices()[*s]; - j next_round_dofs; + + // find all neighbors of the + // dofs numbered in the last + // round + for (unsigned int i=0; i new_entry (coordination, *s); - dofs_by_coordination.insert (new_entry); - }; + next_round_dofs.push_back (sparsity.get_column_numbers()[j]); - //// - std::multimap::iterator i; - for (i = dofs_by_coordination.begin(); i!=dofs_by_coordination.end(); ++i) - new_indices[i->second] = next_free_number++; - - // after that: copy this round's - // dofs for the next round - last_round_dofs = next_round_dofs; - }; + // sort dof numbers + std::sort (next_round_dofs.begin(), next_round_dofs.end()); + + // delete multiple entries + std::vector::iterator end_sorted; + end_sorted = std::unique (next_round_dofs.begin(), next_round_dofs.end()); + next_round_dofs.erase (end_sorted, next_round_dofs.end()); + + // eliminate dofs which are + // already numbered + for (int s=next_round_dofs.size()-1; s>=0; --s) + if (new_indices[next_round_dofs[s]] != DoFHandler::invalid_dof_index) + next_round_dofs.erase (next_round_dofs.begin() + s); + + // check whether there are any new + // dofs + all_dofs_renumbered = (next_round_dofs.size() == 0); + if (all_dofs_renumbered) + // end loop if possible + continue; + + + // store for each coordination + // number the dofs with these + // coordination number + std::multimap dofs_by_coordination; + + // find coordination number for + // each of these dofs + for (std::vector::iterator s=next_round_dofs.begin(); + s!=next_round_dofs.end(); ++s) + { + unsigned int coordination = 0; + for (unsigned int j=sparsity.get_rowstart_indices()[*s]; + j new_entry (coordination, *s); + dofs_by_coordination.insert (new_entry); + }; + + //// + std::multimap::iterator i; + for (i = dofs_by_coordination.begin(); i!=dofs_by_coordination.end(); ++i) + new_indices[i->second] = next_free_number++; + + // after that: copy this round's + // dofs for the next round + last_round_dofs = next_round_dofs; + }; #ifdef DEBUG - // test for all indices numbered + // test for all indices numbered //TODO: Test fails. Do not check before unification. - if (std::find (new_indices.begin(), new_indices.end(), - DoFHandler::invalid_dof_index) - != - new_indices.end()) - Assert (false, ExcRenumberingIncomplete()); - Assert (next_free_number == n_dofs, - ExcRenumberingIncomplete()); + if (std::find (new_indices.begin(), new_indices.end(), + DoFHandler::invalid_dof_index) + != + new_indices.end()) + Assert (false, ExcRenumberingIncomplete()); + Assert (next_free_number == n_dofs, + ExcRenumberingIncomplete()); #endif - if (reversed_numbering) - for (std::vector::iterator i=new_indices.begin(); i!=new_indices.end(); ++i) - *i = n_dofs-*i; + if (reversed_numbering) + for (std::vector::iterator i=new_indices.begin(); i!=new_indices.end(); ++i) + *i = n_dofs-*i; - // actually perform renumbering; - // this is dimension specific and - // thus needs an own function - dof_handler.renumber_dofs (level, new_indices); -} + // actually perform renumbering; + // this is dimension specific and + // thus needs an own function + dof_handler.renumber_dofs (level, new_indices); + } -template -void -DoFRenumbering:: -component_wise (DoFHandler &dof_handler, - const std::vector &component_order_arg) -{ - std::vector renumbering (dof_handler.n_dofs(), - DoFHandler::invalid_dof_index); - - typedef - WrapDoFIterator::active_cell_iterator> ITERATOR; - - typename DoFHandler::active_cell_iterator - istart = dof_handler.begin_active(); - ITERATOR start = istart; - const typename DoFHandler::cell_iterator - end = dof_handler.end(); - - unsigned int result = - compute_component_wise::cell_iterator>(renumbering, - start, end, - component_order_arg); + template + void + + component_wise (DoFHandler &dof_handler, + const std::vector &component_order_arg) + { + std::vector renumbering (dof_handler.n_dofs(), + DoFHandler::invalid_dof_index); + + typedef + internal::WrapDoFIterator::active_cell_iterator> ITERATOR; + + typename DoFHandler::active_cell_iterator + istart = dof_handler.begin_active(); + ITERATOR start = istart; + const typename DoFHandler::cell_iterator + end = dof_handler.end(); - if (result == 0) return; + unsigned int result = + compute_component_wise::cell_iterator>(renumbering, + start, end, + component_order_arg); + + if (result == 0) return; - Assert (result == dof_handler.n_dofs(), - ExcRenumberingIncomplete()); + Assert (result == dof_handler.n_dofs(), + ExcRenumberingIncomplete()); - dof_handler.renumber_dofs (renumbering); -} + dof_handler.renumber_dofs (renumbering); + } -template -void DoFRenumbering:: -component_wise (MGDoFHandler &dof_handler, - const unsigned int level, - const std::vector &component_order_arg) -{ - std::vector renumbering (dof_handler.n_dofs(level), - DoFHandler::invalid_dof_index); + template + void + component_wise (MGDoFHandler &dof_handler, + const unsigned int level, + const std::vector &component_order_arg) + { + std::vector renumbering (dof_handler.n_dofs(level), + DoFHandler::invalid_dof_index); - typedef - WrapMGDoFIterator::cell_iterator> ITERATOR; + typedef + internal::WrapMGDoFIterator::cell_iterator> ITERATOR; - typename MGDoFHandler::cell_iterator - istart =dof_handler.begin(level); - ITERATOR start = istart; - typename MGDoFHandler::cell_iterator - iend = dof_handler.end(level); - const ITERATOR end = iend; + typename MGDoFHandler::cell_iterator + istart =dof_handler.begin(level); + ITERATOR start = istart; + typename MGDoFHandler::cell_iterator + iend = dof_handler.end(level); + const ITERATOR end = iend; - const unsigned int result = - compute_component_wise( - renumbering, start, end, component_order_arg); + const unsigned int result = + compute_component_wise( + renumbering, start, end, component_order_arg); - if (result == 0) return; + if (result == 0) return; - Assert (result == dof_handler.n_dofs(level), - ExcRenumberingIncomplete()); + Assert (result == dof_handler.n_dofs(level), + ExcRenumberingIncomplete()); - if (renumbering.size()!=0) - dof_handler.renumber_dofs (level, renumbering); -} + if (renumbering.size()!=0) + dof_handler.renumber_dofs (level, renumbering); + } -template -void DoFRenumbering:: -component_wise (MGDoFHandler &dof_handler, - const std::vector &component_order_arg) -{ - Threads::ThreadGroup<> threads; + template + void + component_wise (MGDoFHandler &dof_handler, + const std::vector &component_order_arg) + { + Threads::ThreadGroup<> threads; - void (*non_mg_part) (DoFHandler &, const std::vector &) - = &DoFRenumbering::component_wise; - void (*mg_part) (MGDoFHandler &, unsigned int, const std::vector &) - = &DoFRenumbering::component_wise; + void (*non_mg_part) (DoFHandler &, const std::vector &) + = &component_wise; + void (*mg_part) (MGDoFHandler &, unsigned int, const std::vector &) + = &component_wise; - threads += Threads::spawn (non_mg_part) (static_cast&> (dof_handler), - component_order_arg); - for (unsigned int level=0; level&> (dof_handler), + component_order_arg); + for (unsigned int level=0; level -unsigned int -DoFRenumbering:: -compute_component_wise (std::vector& new_indices, - ITERATOR& start, - const ENDITERATOR& end, - const std::vector &component_order_arg) -{ + template + unsigned int + compute_component_wise (std::vector& new_indices, + ITERATOR& start, + const ENDITERATOR& end, + const std::vector &component_order_arg) + { //TODO: Modify for hp - const FiniteElement& fe = start->get_fe(); - const unsigned int dofs_per_cell = fe.dofs_per_cell; + const FiniteElement& fe = start->get_fe(); + const unsigned int dofs_per_cell = fe.dofs_per_cell; // Assert (new_indices.size() == start->dof_handler().n_dofs(), // ExcDimensionMismatch(new_indices.size(), // start->dof_handler.n_dofs())); - // do nothing if the FE has only - // one component - if (fe.n_components() == 1) - { - new_indices.resize(0); - return 0; - } + // do nothing if the FE has only + // one component + if (fe.n_components() == 1) + { + new_indices.resize(0); + return 0; + } - // Copy last argument into a - // writable vector. - std::vector component_order (component_order_arg); - // If the last argument was an - // empty vector, set up things to - // store components in the order - // found in the system. - if (component_order.size() == 0) - for (unsigned int i=0; i local_dof_indices(dofs_per_cell); - - // prebuilt list to which component - // a given dof on a cell - // should go. note that we get into - // trouble here if the shape - // function is not primitive, since - // then there is no single vector - // component to which it - // belongs. in this case, assign it - // to the first vector component to - // which it belongs - std::vector component_list (dofs_per_cell); - for (unsigned int i=0; i component_order (component_order_arg); + // If the last argument was an + // empty vector, set up things to + // store components in the order + // found in the system. + if (component_order.size() == 0) + for (unsigned int i=0; i local_dof_indices(dofs_per_cell); + + // prebuilt list to which component + // a given dof on a cell + // should go. note that we get into + // trouble here if the shape + // function is not primitive, since + // then there is no single vector + // component to which it + // belongs. in this case, assign it + // to the first vector component to + // which it belongs + std::vector component_list (dofs_per_cell); + for (unsigned int i=0; i > + component_to_dof_map (fe.n_components()); + for (;start!=end;++start) { - const unsigned int base_element = - fe.system_to_base_index(i).first.first; - const unsigned int base_multiplicity = - fe.system_to_base_index(i).first.second; - // sum up the number of - // components all the - // elements before that have - unsigned int c=0; - for (unsigned int b=0; b > - component_to_dof_map (fe.n_components()); - for (;start!=end;++start) - { - // on each cell: get dof indices - // and insert them into the global - // list using their component - start.get_dof_indices (local_dof_indices); - for (unsigned int i=0; i::const_iterator - begin_of_component = component_to_dof_map[component].begin(), - end_of_component = component_to_dof_map[component].end(); + const typename std::vector::const_iterator + begin_of_component = component_to_dof_map[component].begin(), + end_of_component = component_to_dof_map[component].end(); - for (typename std::vector::const_iterator - dof_index = begin_of_component; - dof_index != end_of_component; ++dof_index) - new_indices[*dof_index] = next_free_index++; - }; + for (typename std::vector::const_iterator + dof_index = begin_of_component; + dof_index != end_of_component; ++dof_index) + new_indices[*dof_index] = next_free_index++; + }; - return next_free_index; -} + return next_free_index; + } -template -void -DoFRenumbering::sort_selected_dofs_back ( - DH& dof_handler, - const std::vector& selected_dofs) -{ - std::vector renumbering(dof_handler.n_dofs(), - DH::invalid_dof_index); - compute_sort_selected_dofs_back(renumbering, dof_handler, selected_dofs); + template + void + sort_selected_dofs_back (DH& dof_handler, + const std::vector& selected_dofs) + { + std::vector renumbering(dof_handler.n_dofs(), + DH::invalid_dof_index); + compute_sort_selected_dofs_back(renumbering, dof_handler, selected_dofs); - dof_handler.renumber_dofs(renumbering); -} + dof_handler.renumber_dofs(renumbering); + } -template -void -DoFRenumbering::compute_sort_selected_dofs_back ( - std::vector& new_indices, - const DH& dof_handler, - const std::vector& selected_dofs) -{ - const unsigned int n_dofs = dof_handler.n_dofs(); - Assert (selected_dofs.size() == n_dofs, - ExcDimensionMismatch (selected_dofs.size(), n_dofs)); - - // re-sort the dofs according to - // their selection state - Assert (new_indices.size() == n_dofs, - ExcDimensionMismatch(new_indices.size(), n_dofs)); - - const unsigned int n_selected_dofs = std::count (selected_dofs.begin(), - selected_dofs.end(), - false); - - unsigned int next_unselected = 0; - unsigned int next_selected = n_selected_dofs; - for (unsigned int i=0; i + void + compute_sort_selected_dofs_back (std::vector& new_indices, + const DH& dof_handler, + const std::vector& selected_dofs) + { + const unsigned int n_dofs = dof_handler.n_dofs(); + Assert (selected_dofs.size() == n_dofs, + ExcDimensionMismatch (selected_dofs.size(), n_dofs)); + // re-sort the dofs according to + // their selection state + Assert (new_indices.size() == n_dofs, + ExcDimensionMismatch(new_indices.size(), n_dofs)); + + const unsigned int n_selected_dofs = std::count (selected_dofs.begin(), + selected_dofs.end(), + false); + + unsigned int next_unselected = 0; + unsigned int next_selected = n_selected_dofs; + for (unsigned int i=0; i -void -DoFRenumbering::cell_wise_dg ( - DH& dof, - const std::vector& cells) -{ - std::vector renumbering(dof.n_dofs()); - std::vector reverse(dof.n_dofs()); - compute_cell_wise_dg(renumbering, reverse, dof, cells); + + template + void + cell_wise_dg (DH& dof, + const std::vector& cells) + { + std::vector renumbering(dof.n_dofs()); + std::vector reverse(dof.n_dofs()); + compute_cell_wise_dg(renumbering, reverse, dof, cells); - dof.renumber_dofs(renumbering); -} + dof.renumber_dofs(renumbering); + } + -template -void -DoFRenumbering::cell_wise ( - DH& dof, - const std::vector& cells) -{ - std::vector renumbering(dof.n_dofs()); - std::vector reverse(dof.n_dofs()); - compute_cell_wise(renumbering, reverse, dof, cells); + template + void + cell_wise ( + DH& dof, + const std::vector& cells) + { + std::vector renumbering(dof.n_dofs()); + std::vector reverse(dof.n_dofs()); + compute_cell_wise(renumbering, reverse, dof, cells); - dof.renumber_dofs(renumbering); -} + dof.renumber_dofs(renumbering); + } //TODO: Discuss if we cannot replace this function by the next -template -void -DoFRenumbering::compute_cell_wise_dg ( - std::vector& new_indices, - std::vector& reverse, - const DH& dof, - const typename std::vector& cells) -{ - Assert(cells.size() == dof.get_tria().n_active_cells(), - ExcDimensionMismatch(cells.size(), - dof.get_tria().n_active_cells())); + template + void + compute_cell_wise_dg ( + std::vector& new_indices, + std::vector& reverse, + const DH& dof, + const typename std::vector& cells) + { + Assert(cells.size() == dof.get_tria().n_active_cells(), + ExcDimensionMismatch(cells.size(), + dof.get_tria().n_active_cells())); - unsigned int n_global_dofs = dof.n_dofs(); + unsigned int n_global_dofs = dof.n_dofs(); - // Actually, we compute the inverse of the reordering vector, called reverse here. - // Later, irs inverse is computed into new_indices, which is the return argument. + // Actually, we compute the + // inverse of the reordering + // vector, called reverse here. + // Later, irs inverse is computed + // into new_indices, which is the + // return argument. - Assert(new_indices.size() == n_global_dofs, - ExcDimensionMismatch(new_indices.size(), n_global_dofs)); - Assert(reverse.size() == n_global_dofs, - ExcDimensionMismatch(reverse.size(), n_global_dofs)); + Assert(new_indices.size() == n_global_dofs, + ExcDimensionMismatch(new_indices.size(), n_global_dofs)); + Assert(reverse.size() == n_global_dofs, + ExcDimensionMismatch(reverse.size(), n_global_dofs)); - std::vector cell_dofs; + std::vector cell_dofs; - unsigned int global_index = 0; + unsigned int global_index = 0; - typename std::vector::const_iterator cell; + typename std::vector::const_iterator cell; - for(cell = cells.begin(); cell != cells.end(); ++cell) - { - Assert((*cell)->get_fe().n_dofs_per_face()==0, ExcNotDGFEM()); - // Determine the number of dofs - // on this cell and reinit the - // vector storing these - // numbers. - unsigned int n_cell_dofs = (*cell)->get_fe().n_dofs_per_cell(); - cell_dofs.resize(n_cell_dofs); + for(cell = cells.begin(); cell != cells.end(); ++cell) + { + Assert((*cell)->get_fe().n_dofs_per_face()==0, ExcNotDGFEM()); + // Determine the number of dofs + // on this cell and reinit the + // vector storing these + // numbers. + unsigned int n_cell_dofs = (*cell)->get_fe().n_dofs_per_cell(); + cell_dofs.resize(n_cell_dofs); - (*cell)->get_dof_indices(cell_dofs); - - // Sort here to make sure that - // degrees of freedom inside a - // single cell are in the same - // order after renumbering. - std::sort(cell_dofs.begin(), cell_dofs.end()); - - for (unsigned int i=0;iget_dof_indices(cell_dofs); + + // Sort here to make sure that + // degrees of freedom inside a + // single cell are in the same + // order after renumbering. + std::sort(cell_dofs.begin(), cell_dofs.end()); + + for (unsigned int i=0;i -void -DoFRenumbering::compute_cell_wise ( - std::vector& new_indices, - std::vector& reverse, - const DH& dof, - const typename std::vector& cells) -{ - Assert(cells.size() == dof.get_tria().n_active_cells(), - ExcDimensionMismatch(cells.size(), - dof.get_tria().n_active_cells())); + template + void + compute_cell_wise ( + std::vector& new_indices, + std::vector& reverse, + const DH& dof, + const typename std::vector& cells) + { + Assert(cells.size() == dof.get_tria().n_active_cells(), + ExcDimensionMismatch(cells.size(), + dof.get_tria().n_active_cells())); - unsigned int n_global_dofs = dof.n_dofs(); + unsigned int n_global_dofs = dof.n_dofs(); - // Actually, we compute the inverse of the reordering vector, called reverse here. - // Later, irs inverse is computed into new_indices, which is the return argument. + // Actually, we compute the + // inverse of the reordering + // vector, called reverse here. + // Later, irs inverse is computed + // into new_indices, which is the + // return argument. - Assert(new_indices.size() == n_global_dofs, - ExcDimensionMismatch(new_indices.size(), n_global_dofs)); - Assert(reverse.size() == n_global_dofs, - ExcDimensionMismatch(reverse.size(), n_global_dofs)); + Assert(new_indices.size() == n_global_dofs, + ExcDimensionMismatch(new_indices.size(), n_global_dofs)); + Assert(reverse.size() == n_global_dofs, + ExcDimensionMismatch(reverse.size(), n_global_dofs)); - // For continuous elements, we must - // make sure, that each dof is - // reordered only once. - std::vector already_sorted(n_global_dofs, false); - std::vector cell_dofs; + // For continuous elements, we must + // make sure, that each dof is + // reordered only once. + std::vector already_sorted(n_global_dofs, false); + std::vector cell_dofs; - unsigned int global_index = 0; + unsigned int global_index = 0; - typename std::vector::const_iterator cell; + typename std::vector::const_iterator cell; - for(cell = cells.begin(); cell != cells.end(); ++cell) - { - // Determine the number of dofs - // on this cell and reinit the - // vector storing these - // numbers. - unsigned int n_cell_dofs = (*cell)->get_fe().n_dofs_per_cell(); - cell_dofs.resize(n_cell_dofs); + for(cell = cells.begin(); cell != cells.end(); ++cell) + { + // Determine the number of dofs + // on this cell and reinit the + // vector storing these + // numbers. + unsigned int n_cell_dofs = (*cell)->get_fe().n_dofs_per_cell(); + cell_dofs.resize(n_cell_dofs); - (*cell)->get_dof_indices(cell_dofs); - - // Sort here to make sure that - // degrees of freedom inside a - // single cell are in the same - // order after renumbering. - std::sort(cell_dofs.begin(), cell_dofs.end()); - - for (unsigned int i=0;iget_dof_indices(cell_dofs); + + // Sort here to make sure that + // degrees of freedom inside a + // single cell are in the same + // order after renumbering. + std::sort(cell_dofs.begin(), cell_dofs.end()); + + for (unsigned int i=0;i -void DoFRenumbering::cell_wise_dg ( - MGDoFHandler& dof, - const unsigned int level, - const typename std::vector::cell_iterator>& cells) -{ - std::vector renumbering(dof.n_dofs(level)); - std::vector reverse(dof.n_dofs(level)); + template + void cell_wise_dg ( + MGDoFHandler& dof, + const unsigned int level, + const typename std::vector::cell_iterator>& cells) + { + std::vector renumbering(dof.n_dofs(level)); + std::vector reverse(dof.n_dofs(level)); - compute_cell_wise_dg(renumbering, reverse, dof, level, cells); - dof.renumber_dofs(level, renumbering); -} + compute_cell_wise_dg(renumbering, reverse, dof, level, cells); + dof.renumber_dofs(level, renumbering); + } + -template -void DoFRenumbering::cell_wise ( - MGDoFHandler& dof, - const unsigned int level, - const typename std::vector::cell_iterator>& cells) -{ - std::vector renumbering(dof.n_dofs(level)); - std::vector reverse(dof.n_dofs(level)); + template + void cell_wise ( + MGDoFHandler& dof, + const unsigned int level, + const typename std::vector::cell_iterator>& cells) + { + std::vector renumbering(dof.n_dofs(level)); + std::vector reverse(dof.n_dofs(level)); - compute_cell_wise(renumbering, reverse, dof, level, cells); - dof.renumber_dofs(level, renumbering); -} + compute_cell_wise(renumbering, reverse, dof, level, cells); + dof.renumber_dofs(level, renumbering); + } + -template -void DoFRenumbering::compute_cell_wise_dg ( - std::vector& new_order, - std::vector& reverse, - const MGDoFHandler& dof, - const unsigned int level, - const typename std::vector::cell_iterator>& cells) -{ - Assert(cells.size() == dof.get_tria().n_cells(level), - ExcDimensionMismatch(cells.size(), - dof.get_tria().n_cells(level))); - switch (dim) - { - case 3: - Assert(dof.get_fe().n_dofs_per_quad()==0, - ExcNotDGFEM()); - case 2: - Assert(dof.get_fe().n_dofs_per_line()==0, - ExcNotDGFEM()); - default: - Assert(dof.get_fe().n_dofs_per_vertex()==0, - ExcNotDGFEM()); - } + template + void compute_cell_wise_dg ( + std::vector& new_order, + std::vector& reverse, + const MGDoFHandler& dof, + const unsigned int level, + const typename std::vector::cell_iterator>& cells) + { + Assert(cells.size() == dof.get_tria().n_cells(level), + ExcDimensionMismatch(cells.size(), + dof.get_tria().n_cells(level))); + switch (dim) + { + case 3: + Assert(dof.get_fe().n_dofs_per_quad()==0, + ExcNotDGFEM()); + case 2: + Assert(dof.get_fe().n_dofs_per_line()==0, + ExcNotDGFEM()); + default: + Assert(dof.get_fe().n_dofs_per_vertex()==0, + ExcNotDGFEM()); + } - Assert (new_order.size() == dof.n_dofs(level), - ExcDimensionMismatch(new_order.size(), dof.n_dofs(level))); - Assert (reverse.size() == dof.n_dofs(level), - ExcDimensionMismatch(reverse.size(), dof.n_dofs(level))); + Assert (new_order.size() == dof.n_dofs(level), + ExcDimensionMismatch(new_order.size(), dof.n_dofs(level))); + Assert (reverse.size() == dof.n_dofs(level), + ExcDimensionMismatch(reverse.size(), dof.n_dofs(level))); - unsigned int n_global_dofs = dof.n_dofs(level); - unsigned int n_cell_dofs = dof.get_fe().n_dofs_per_cell(); + unsigned int n_global_dofs = dof.n_dofs(level); + unsigned int n_cell_dofs = dof.get_fe().n_dofs_per_cell(); - std::vector cell_dofs(n_cell_dofs); + std::vector cell_dofs(n_cell_dofs); - unsigned int global_index = 0; + unsigned int global_index = 0; - typename std::vector::cell_iterator>::const_iterator cell; + typename std::vector::cell_iterator>::const_iterator cell; - for(cell = cells.begin(); cell != cells.end(); ++cell) - { - Assert ((*cell)->level() == (int) level, ExcInternalError()); + for(cell = cells.begin(); cell != cells.end(); ++cell) + { + Assert ((*cell)->level() == (int) level, ExcInternalError()); - (*cell)->get_mg_dof_indices(cell_dofs); - std::sort(cell_dofs.begin(), cell_dofs.end()); + (*cell)->get_mg_dof_indices(cell_dofs); + std::sort(cell_dofs.begin(), cell_dofs.end()); - for (unsigned int i=0;i -void DoFRenumbering::compute_cell_wise ( - std::vector& new_order, - std::vector& reverse, - const MGDoFHandler& dof, - const unsigned int level, - const typename std::vector::cell_iterator>& cells) -{ - Assert(cells.size() == dof.get_tria().n_cells(level), - ExcDimensionMismatch(cells.size(), - dof.get_tria().n_cells(level))); - Assert (new_order.size() == dof.n_dofs(level), - ExcDimensionMismatch(new_order.size(), dof.n_dofs(level))); - Assert (reverse.size() == dof.n_dofs(level), - ExcDimensionMismatch(reverse.size(), dof.n_dofs(level))); + template + void compute_cell_wise ( + std::vector& new_order, + std::vector& reverse, + const MGDoFHandler& dof, + const unsigned int level, + const typename std::vector::cell_iterator>& cells) + { + Assert(cells.size() == dof.get_tria().n_cells(level), + ExcDimensionMismatch(cells.size(), + dof.get_tria().n_cells(level))); + Assert (new_order.size() == dof.n_dofs(level), + ExcDimensionMismatch(new_order.size(), dof.n_dofs(level))); + Assert (reverse.size() == dof.n_dofs(level), + ExcDimensionMismatch(reverse.size(), dof.n_dofs(level))); - unsigned int n_global_dofs = dof.n_dofs(level); - unsigned int n_cell_dofs = dof.get_fe().n_dofs_per_cell(); + unsigned int n_global_dofs = dof.n_dofs(level); + unsigned int n_cell_dofs = dof.get_fe().n_dofs_per_cell(); - std::vector already_sorted(n_global_dofs, false); - std::vector cell_dofs(n_cell_dofs); + std::vector already_sorted(n_global_dofs, false); + std::vector cell_dofs(n_cell_dofs); - unsigned int global_index = 0; + unsigned int global_index = 0; - typename std::vector::cell_iterator>::const_iterator cell; + typename std::vector::cell_iterator>::const_iterator cell; - for(cell = cells.begin(); cell != cells.end(); ++cell) - { - Assert ((*cell)->level() == (int) level, ExcInternalError()); + for(cell = cells.begin(); cell != cells.end(); ++cell) + { + Assert ((*cell)->level() == (int) level, ExcInternalError()); - (*cell)->get_mg_dof_indices(cell_dofs); - std::sort(cell_dofs.begin(), cell_dofs.end()); + (*cell)->get_mg_dof_indices(cell_dofs); + std::sort(cell_dofs.begin(), cell_dofs.end()); - for (unsigned int i=0;i -void -DoFRenumbering::downstream_dg (DH& dof, const Point& direction) -{ - std::vector renumbering(dof.n_dofs()); - compute_downstream_dg(renumbering, dof, direction); + template + void + downstream_dg (DH& dof, const Point& direction) + { + std::vector renumbering(dof.n_dofs()); + compute_downstream_dg(renumbering, dof, direction); - dof.renumber_dofs(renumbering); -} + dof.renumber_dofs(renumbering); + } -template -void -DoFRenumbering::downstream (DH& dof, const Point& direction) -{ - std::vector renumbering(dof.n_dofs()); - std::vector reverse(dof.n_dofs()); - compute_downstream(renumbering, reverse, dof, direction); + template + void + downstream (DH& dof, const Point& direction) + { + std::vector renumbering(dof.n_dofs()); + std::vector reverse(dof.n_dofs()); + compute_downstream(renumbering, reverse, dof, direction); - dof.renumber_dofs(renumbering); -} + dof.renumber_dofs(renumbering); + } -template -void -DoFRenumbering::compute_downstream_dg ( - std::vector& new_indices, - const DH& dof, - const Point& direction) -{ - std::vector - ordered_cells(dof.get_tria().n_active_cells()); - const CompareDownstream comparator(direction); + template + void + compute_downstream_dg ( + std::vector& new_indices, + const DH& dof, + const Point& direction) + { + std::vector + ordered_cells(dof.get_tria().n_active_cells()); + const internal::CompareDownstream comparator(direction); - typename DH::active_cell_iterator begin = dof.begin_active(); - typename DH::active_cell_iterator end = dof.end(); + typename DH::active_cell_iterator begin = dof.begin_active(); + typename DH::active_cell_iterator end = dof.end(); - copy (begin, end, ordered_cells.begin()); - std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); + copy (begin, end, ordered_cells.begin()); + std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); - std::vector reverse(new_indices.size()); - compute_cell_wise_dg(new_indices, reverse, dof, ordered_cells); -} + std::vector reverse(new_indices.size()); + compute_cell_wise_dg(new_indices, reverse, dof, ordered_cells); + } + -template -void -DoFRenumbering::compute_downstream_dg ( - std::vector& new_indices, - std::vector& reverse, - const DH& dof, - const Point& direction) -{ - std::vector - ordered_cells(dof.get_tria().n_active_cells()); - const CompareDownstream comparator(direction); + template + void + compute_downstream_dg ( + std::vector& new_indices, + std::vector& reverse, + const DH& dof, + const Point& direction) + { + std::vector + ordered_cells(dof.get_tria().n_active_cells()); + const internal::CompareDownstream comparator(direction); - typename DH::active_cell_iterator begin = dof.begin_active(); - typename DH::active_cell_iterator end = dof.end(); + typename DH::active_cell_iterator begin = dof.begin_active(); + typename DH::active_cell_iterator end = dof.end(); - copy (begin, end, ordered_cells.begin()); - std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); + copy (begin, end, ordered_cells.begin()); + std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); - compute_cell_wise_dg(new_indices, reverse, dof, ordered_cells); -} + compute_cell_wise_dg(new_indices, reverse, dof, ordered_cells); + } -template -void -DoFRenumbering::compute_downstream ( - std::vector& new_indices, - std::vector& reverse, - const DH& dof, - const Point& direction) -{ - std::vector - ordered_cells(dof.get_tria().n_active_cells()); - const CompareDownstream comparator(direction); + template + void + compute_downstream ( + std::vector& new_indices, + std::vector& reverse, + const DH& dof, + const Point& direction) + { + std::vector + ordered_cells(dof.get_tria().n_active_cells()); + const internal::CompareDownstream comparator(direction); - typename DH::active_cell_iterator begin = dof.begin_active(); - typename DH::active_cell_iterator end = dof.end(); + typename DH::active_cell_iterator begin = dof.begin_active(); + typename DH::active_cell_iterator end = dof.end(); - copy (begin, end, ordered_cells.begin()); - std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); + copy (begin, end, ordered_cells.begin()); + std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); - compute_cell_wise(new_indices, reverse, dof, ordered_cells); -} + compute_cell_wise(new_indices, reverse, dof, ordered_cells); + } -template -void DoFRenumbering::downstream_dg (MGDoFHandler& dof, - const unsigned int level, - const Point& direction) -{ - std::vector renumbering(dof.n_dofs(level)); - std::vector reverse(dof.n_dofs(level)); - compute_downstream_dg(renumbering, reverse, dof, level, direction); - dof.renumber_dofs(level, renumbering); -} + template + void downstream_dg (MGDoFHandler& dof, + const unsigned int level, + const Point& direction) + { + std::vector renumbering(dof.n_dofs(level)); + std::vector reverse(dof.n_dofs(level)); + compute_downstream_dg(renumbering, reverse, dof, level, direction); + + dof.renumber_dofs(level, renumbering); + } + -template -void DoFRenumbering::downstream (MGDoFHandler& dof, - const unsigned int level, - const Point& direction) -{ - std::vector renumbering(dof.n_dofs(level)); - std::vector reverse(dof.n_dofs(level)); - compute_downstream(renumbering, reverse, dof, level, direction); + template + void downstream (MGDoFHandler& dof, + const unsigned int level, + const Point& direction) + { + std::vector renumbering(dof.n_dofs(level)); + std::vector reverse(dof.n_dofs(level)); + compute_downstream(renumbering, reverse, dof, level, direction); - dof.renumber_dofs(level, renumbering); -} + dof.renumber_dofs(level, renumbering); + } + -template -void -DoFRenumbering::compute_downstream_dg ( - std::vector& new_indices, - std::vector& reverse, - const MGDoFHandler& dof, - const unsigned int level, - const Point& direction) -{ - std::vector::cell_iterator> - ordered_cells(dof.get_tria().n_cells(level)); - const CompareDownstream comparator(direction); + template + void + compute_downstream_dg ( + std::vector& new_indices, + std::vector& reverse, + const MGDoFHandler& dof, + const unsigned int level, + const Point& direction) + { + std::vector::cell_iterator> + ordered_cells(dof.get_tria().n_cells(level)); + const internal::CompareDownstream comparator(direction); - typename MGDoFHandler::cell_iterator begin = dof.begin(level); - typename MGDoFHandler::cell_iterator end = dof.end(level); + typename MGDoFHandler::cell_iterator begin = dof.begin(level); + typename MGDoFHandler::cell_iterator end = dof.end(level); - std::copy (begin, end, ordered_cells.begin()); - std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); + std::copy (begin, end, ordered_cells.begin()); + std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); - compute_cell_wise_dg(new_indices, reverse, dof, level, ordered_cells); -} + compute_cell_wise_dg(new_indices, reverse, dof, level, ordered_cells); + } + -template -void -DoFRenumbering::compute_downstream ( - std::vector& new_indices, - std::vector& reverse, - const MGDoFHandler& dof, - const unsigned int level, - const Point& direction) -{ - std::vector::cell_iterator> - ordered_cells(dof.get_tria().n_cells(level)); - const CompareDownstream comparator(direction); + template + void + compute_downstream ( + std::vector& new_indices, + std::vector& reverse, + const MGDoFHandler& dof, + const unsigned int level, + const Point& direction) + { + std::vector::cell_iterator> + ordered_cells(dof.get_tria().n_cells(level)); + const internal::CompareDownstream comparator(direction); - typename MGDoFHandler::cell_iterator begin = dof.begin(level); - typename MGDoFHandler::cell_iterator end = dof.end(level); + typename MGDoFHandler::cell_iterator begin = dof.begin(level); + typename MGDoFHandler::cell_iterator end = dof.end(level); - std::copy (begin, end, ordered_cells.begin()); - std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); + std::copy (begin, end, ordered_cells.begin()); + std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); - compute_cell_wise(new_indices, reverse, dof, level, ordered_cells); -} + compute_cell_wise(new_indices, reverse, dof, level, ordered_cells); + } /** * Provide comparator for DoFCellAccessors */ - -template -struct ClockCells -{ - /** - * Center of rotation. - */ - const Point& center; - /** - * Revert sorting order. - */ - bool counter; + namespace internal + { + template + struct ClockCells + { + /** + * Center of rotation. + */ + const Point& center; + /** + * Revert sorting order. + */ + bool counter; - /** - * Constructor. - */ - ClockCells (const Point& center, bool counter) : - center(center), - counter(counter) - {} - /** - * Comparison operator - */ - template - bool operator () (const DHCellIterator& c1, - const DHCellIterator& c2) const - { + /** + * Constructor. + */ + ClockCells (const Point& center, bool counter) : + center(center), + counter(counter) + {} + /** + * Comparison operator + */ + template + bool operator () (const DHCellIterator& c1, + const DHCellIterator& c2) const + { - const Point v1 = c1->center() - center; - const Point v2 = c2->center() - center; - const double s1 = std::atan2(v1(0), v1(1)); - const double s2 = std::atan2(v2(0), v2(1)); - return ( counter ? (s1>s2) : (s2>s1)); - } -}; + const Point v1 = c1->center() - center; + const Point v2 = c2->center() - center; + const double s1 = std::atan2(v1(0), v1(1)); + const double s2 = std::atan2(v2(0), v2(1)); + return ( counter ? (s1>s2) : (s2>s1)); + } + }; + } + -template -void -DoFRenumbering::clockwise_dg ( - DH& dof, - const Point& center, - const bool counter) -{ - std::vector renumbering(dof.n_dofs()); - compute_clockwise_dg(renumbering, dof, center, counter); + template + void + clockwise_dg ( + DH& dof, + const Point& center, + const bool counter) + { + std::vector renumbering(dof.n_dofs()); + compute_clockwise_dg(renumbering, dof, center, counter); - dof.renumber_dofs(renumbering); -} + dof.renumber_dofs(renumbering); + } -template -void -DoFRenumbering::compute_clockwise_dg ( - std::vector& new_indices, - const DH& dof, - const Point& center, - const bool counter) -{ - std::vector - ordered_cells(dof.get_tria().n_active_cells()); - ClockCells comparator(center, counter); + template + void + compute_clockwise_dg ( + std::vector& new_indices, + const DH& dof, + const Point& center, + const bool counter) + { + std::vector + ordered_cells(dof.get_tria().n_active_cells()); + internal::ClockCells comparator(center, counter); - typename DH::active_cell_iterator begin = dof.begin_active(); - typename DH::active_cell_iterator end = dof.end(); + typename DH::active_cell_iterator begin = dof.begin_active(); + typename DH::active_cell_iterator end = dof.end(); - std::copy (begin, end, ordered_cells.begin()); - std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); + std::copy (begin, end, ordered_cells.begin()); + std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); - std::vector reverse(new_indices.size()); - compute_cell_wise_dg(new_indices, reverse, dof, ordered_cells); -} + std::vector reverse(new_indices.size()); + compute_cell_wise_dg(new_indices, reverse, dof, ordered_cells); + } + -template -void DoFRenumbering::clockwise_dg (MGDoFHandler& dof, - const unsigned int level, - const Point& center, - const bool counter) -{ - std::vector::cell_iterator> - ordered_cells(dof.get_tria().n_cells(level)); - ClockCells comparator(center, counter); + template + void clockwise_dg (MGDoFHandler& dof, + const unsigned int level, + const Point& center, + const bool counter) + { + std::vector::cell_iterator> + ordered_cells(dof.get_tria().n_cells(level)); + internal::ClockCells comparator(center, counter); - typename MGDoFHandler::cell_iterator begin = dof.begin(level); - typename MGDoFHandler::cell_iterator end = dof.end(level); + typename MGDoFHandler::cell_iterator begin = dof.begin(level); + typename MGDoFHandler::cell_iterator end = dof.end(level); - std::copy (begin, end, ordered_cells.begin()); - std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); + std::copy (begin, end, ordered_cells.begin()); + std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); - cell_wise_dg(dof, level, ordered_cells); -} + cell_wise_dg(dof, level, ordered_cells); + } -template -void -DoFRenumbering::random (DH& dof_handler) -{ - std::vector renumbering(dof_handler.n_dofs(), - DH::invalid_dof_index); - compute_random(renumbering, dof_handler); + template + void + random (DH& dof_handler) + { + std::vector renumbering(dof_handler.n_dofs(), + DH::invalid_dof_index); + compute_random(renumbering, dof_handler); - dof_handler.renumber_dofs(renumbering); -} + dof_handler.renumber_dofs(renumbering); + } + -template -void -DoFRenumbering::compute_random ( - std::vector& new_indices, - const DH& dof_handler) -{ - const unsigned int n_dofs = dof_handler.n_dofs(); - Assert(new_indices.size() == n_dofs, - ExcDimensionMismatch(new_indices.size(), n_dofs)); + template + void + compute_random ( + std::vector& new_indices, + const DH& dof_handler) + { + const unsigned int n_dofs = dof_handler.n_dofs(); + Assert(new_indices.size() == n_dofs, + ExcDimensionMismatch(new_indices.size(), n_dofs)); - for (unsigned i=0; i -void -DoFRenumbering::subdomain_wise (DoFHandler &dof_handler) -{ - std::vector renumbering(dof_handler.n_dofs(), - DoFHandler::invalid_dof_index); - compute_subdomain_wise(renumbering, dof_handler); + template + void + subdomain_wise (DoFHandler &dof_handler) + { + std::vector renumbering(dof_handler.n_dofs(), + DoFHandler::invalid_dof_index); + compute_subdomain_wise(renumbering, dof_handler); - dof_handler.renumber_dofs(renumbering); -} + dof_handler.renumber_dofs(renumbering); + } + -template -void -DoFRenumbering:: -compute_subdomain_wise (std::vector &new_dof_indices, - const DoFHandler &dof_handler) -{ - const unsigned int n_dofs = dof_handler.n_dofs(); - Assert (new_dof_indices.size() == n_dofs, - ExcDimensionMismatch (new_dof_indices.size(), n_dofs)); - - // first get the association of each dof - // with a subdomain and determine the total - // number of subdomain ids used - std::vector subdomain_association (n_dofs); - DoFTools::get_subdomain_association (dof_handler, - subdomain_association); - const unsigned int n_subdomains - = *std::max_element (subdomain_association.begin(), - subdomain_association.end()) + 1; - - // then renumber the subdomains by first - // looking at those belonging to subdomain - // 0, then those of subdomain 1, etc. note - // that the algorithm is stable, i.e. if - // two dofs i,j have i + void + compute_subdomain_wise (std::vector &new_dof_indices, + const DoFHandler &dof_handler) + { + const unsigned int n_dofs = dof_handler.n_dofs(); + Assert (new_dof_indices.size() == n_dofs, + ExcDimensionMismatch (new_dof_indices.size(), n_dofs)); + + // first get the association of each dof + // with a subdomain and determine the total + // number of subdomain ids used + std::vector subdomain_association (n_dofs); + DoFTools::get_subdomain_association (dof_handler, + subdomain_association); + const unsigned int n_subdomains + = *std::max_element (subdomain_association.begin(), + subdomain_association.end()) + 1; + + // then renumber the subdomains by first + // looking at those belonging to subdomain + // 0, then those of subdomain 1, etc. note + // that the algorithm is stable, i.e. if + // two dofs i,j have i &, bool, bool); - -template -void -DoFRenumbering::boost::compute_Cuthill_McKee (std::vector &, - const DoFHandler &, bool, bool); - -template -void -DoFRenumbering::boost::king_ordering (DoFHandler &, bool, bool); - -template -void -DoFRenumbering::boost::compute_king_ordering (std::vector &, - const DoFHandler &, bool, bool); - -template -void -DoFRenumbering::boost::minimum_degree (DoFHandler &, bool, bool); - -template -void -DoFRenumbering::boost::compute_minimum_degree (std::vector &, - const DoFHandler &, bool, bool); - + namespace boost + { +//TODO[WB]: also implement the following boost for hp DoFHandlers etc. + template + void + Cuthill_McKee (DoFHandler &, bool, bool); + + template + void + compute_Cuthill_McKee (std::vector &, + const DoFHandler &, bool, bool); + + template + void + king_ordering (DoFHandler &, bool, bool); + + template + void + compute_king_ordering (std::vector &, + const DoFHandler &, bool, bool); + + template + void + minimum_degree (DoFHandler &, bool, bool); + + template + void + compute_minimum_degree (std::vector &, + const DoFHandler &, bool, bool); + } + // non-boost functions: -template -void DoFRenumbering::Cuthill_McKee > -(DoFHandler&, - const bool, - const bool, - const std::vector&); - -template -void -DoFRenumbering::compute_Cuthill_McKee > -(std::vector&, - const DoFHandler&, - const bool, - const bool, - const std::vector&); - - -template -void DoFRenumbering::Cuthill_McKee > -(hp::DoFHandler&, - const bool, - const bool, - const std::vector&); - -template -void -DoFRenumbering::compute_Cuthill_McKee > -(std::vector&, - const hp::DoFHandler&, - const bool, - const bool, - const std::vector&); - - -template -void DoFRenumbering::component_wise -(DoFHandler&, - const std::vector&); + template + void Cuthill_McKee > + (DoFHandler&, + const bool, + const bool, + const std::vector&); + + template + void + compute_Cuthill_McKee > + (std::vector&, + const DoFHandler&, + const bool, + const bool, + const std::vector&); + + + template + void Cuthill_McKee > + (hp::DoFHandler&, + const bool, + const bool, + const std::vector&); + + template + void + compute_Cuthill_McKee > + (std::vector&, + const hp::DoFHandler&, + const bool, + const bool, + const std::vector&); + + + template + void component_wise + (DoFHandler&, + const std::vector&); // template // void -// DoFRenumbering::compute_component_wise +// compute_component_wise // (std::vector&, // const DoFHandler&, // const std::vector&); -template -void DoFRenumbering::component_wise -(MGDoFHandler&, - unsigned int, - const std::vector&); + template + void component_wise + (MGDoFHandler&, + unsigned int, + const std::vector&); -template -void DoFRenumbering::component_wise -(MGDoFHandler&, - const std::vector&); + template + void component_wise + (MGDoFHandler&, + const std::vector&); // DG renumbering for DoFHandler -template -void -DoFRenumbering::cell_wise_dg > -(DoFHandler&, - const std::vector::cell_iterator>&); - -template -void -DoFRenumbering::compute_cell_wise_dg > -(std::vector&, std::vector&, - const DoFHandler&, - const std::vector::cell_iterator>&); - -template -void -DoFRenumbering::downstream_dg > -(DoFHandler&, - const Point&); - -template -void -DoFRenumbering::compute_downstream_dg > -(std::vector&,std::vector&, - const DoFHandler&, - const Point&); -template -void -DoFRenumbering::compute_downstream_dg > -(std::vector&, - const DoFHandler&, - const Point&); - -template void -DoFRenumbering::cell_wise > -(DoFHandler&, - const std::vector::cell_iterator>&); - -template void -DoFRenumbering::compute_cell_wise > -(std::vector&, std::vector&, - const DoFHandler&, - const std::vector::cell_iterator>&); - -template void -DoFRenumbering::downstream > -(DoFHandler&, const Point&); - -template void -DoFRenumbering::compute_downstream > -(std::vector&,std::vector&, - const DoFHandler&, const Point&); - -template -void -DoFRenumbering::clockwise_dg > -(DoFHandler&, const Point&, bool); - -template -void -DoFRenumbering::compute_clockwise_dg > -(std::vector&, const DoFHandler&, - const Point&, const bool); + template + void + cell_wise_dg > + (DoFHandler&, + const std::vector::cell_iterator>&); + + template + void + compute_cell_wise_dg > + (std::vector&, std::vector&, + const DoFHandler&, + const std::vector::cell_iterator>&); + + template + void + downstream_dg > + (DoFHandler&, + const Point&); + + template + void + compute_downstream_dg > + (std::vector&,std::vector&, + const DoFHandler&, + const Point&); + template + void + compute_downstream_dg > + (std::vector&, + const DoFHandler&, + const Point&); + + template void + cell_wise > + (DoFHandler&, + const std::vector::cell_iterator>&); + + template void + compute_cell_wise > + (std::vector&, std::vector&, + const DoFHandler&, + const std::vector::cell_iterator>&); + + template void + downstream > + (DoFHandler&, const Point&); + + template void + compute_downstream > + (std::vector&,std::vector&, + const DoFHandler&, const Point&); + + template + void + clockwise_dg > + (DoFHandler&, const Point&, bool); + + template + void + compute_clockwise_dg > + (std::vector&, const DoFHandler&, + const Point&, const bool); // Renumbering for hp::DoFHandler -template void -DoFRenumbering::cell_wise_dg > -(hp::DoFHandler&, - const std::vector::cell_iterator>&); - -template void -DoFRenumbering::compute_cell_wise_dg > -(std::vector&, std::vector&, - const hp::DoFHandler&, - const std::vector::cell_iterator>&); - -template void -DoFRenumbering::cell_wise > -(hp::DoFHandler&, - const std::vector::cell_iterator>&); - -template void -DoFRenumbering::compute_cell_wise > -(std::vector&, std::vector&, - const hp::DoFHandler&, - const std::vector::cell_iterator>&); - -template void -DoFRenumbering::downstream_dg > -(hp::DoFHandler&, - const Point&); - -template void -DoFRenumbering::compute_downstream_dg > -(std::vector&,std::vector&, - const hp::DoFHandler&, - const Point&); -template -void -DoFRenumbering::compute_downstream_dg > -(std::vector&, - const hp::DoFHandler&, - const Point&); - -template void -DoFRenumbering::downstream > -(hp::DoFHandler&, - const Point&); - -template void -DoFRenumbering::compute_downstream > -(std::vector&,std::vector&, - const hp::DoFHandler&, - const Point&); - -template -void -DoFRenumbering::clockwise_dg > -(hp::DoFHandler&, - const Point&, bool); - -template -void -DoFRenumbering::compute_clockwise_dg > -(std::vector&, - const hp::DoFHandler&, - const Point&, - const bool); + template void + cell_wise_dg > + (hp::DoFHandler&, + const std::vector::cell_iterator>&); + + template void + compute_cell_wise_dg > + (std::vector&, std::vector&, + const hp::DoFHandler&, + const std::vector::cell_iterator>&); + + template void + cell_wise > + (hp::DoFHandler&, + const std::vector::cell_iterator>&); + + template void + compute_cell_wise > + (std::vector&, std::vector&, + const hp::DoFHandler&, + const std::vector::cell_iterator>&); + + template void + downstream_dg > + (hp::DoFHandler&, + const Point&); + + template void + compute_downstream_dg > + (std::vector&,std::vector&, + const hp::DoFHandler&, + const Point&); + template + void + compute_downstream_dg > + (std::vector&, + const hp::DoFHandler&, + const Point&); + + template void + downstream > + (hp::DoFHandler&, + const Point&); + + template void + compute_downstream > + (std::vector&,std::vector&, + const hp::DoFHandler&, + const Point&); + + template + void + clockwise_dg > + (hp::DoFHandler&, + const Point&, bool); + + template + void + compute_clockwise_dg > + (std::vector&, + const hp::DoFHandler&, + const Point&, + const bool); // MG -template -void DoFRenumbering::downstream_dg -(MGDoFHandler&, const unsigned int, - const Point&); + template + void downstream_dg + (MGDoFHandler&, const unsigned int, + const Point&); -template -void DoFRenumbering::downstream -(MGDoFHandler&, const unsigned int, - const Point&); + template + void downstream + (MGDoFHandler&, const unsigned int, + const Point&); -template -void DoFRenumbering::clockwise_dg -(MGDoFHandler&, - const unsigned int, - const Point&, bool); + template + void clockwise_dg + (MGDoFHandler&, + const unsigned int, + const Point&, bool); // Generic numbering schemes -template -void DoFRenumbering::random > -(DoFHandler &); - -template -void -DoFRenumbering::compute_random > -(std::vector&, - const DoFHandler &); - -template -void DoFRenumbering::sort_selected_dofs_back > -(DoFHandler &, - const std::vector &); - -template -void -DoFRenumbering::compute_sort_selected_dofs_back > -(std::vector&, - const DoFHandler &, - const std::vector &); - -template -void DoFRenumbering::random > -(hp::DoFHandler &); - -template -void -DoFRenumbering::compute_random > -(std::vector&, - const hp::DoFHandler &); - -template -void DoFRenumbering::sort_selected_dofs_back > -(hp::DoFHandler &, - const std::vector &); - -template -void -DoFRenumbering::compute_sort_selected_dofs_back > -(std::vector&, - const hp::DoFHandler &, - const std::vector &); - -template -void DoFRenumbering::subdomain_wise -(DoFHandler &); - -template -void -DoFRenumbering::compute_subdomain_wise -(std::vector&, - const DoFHandler &); - -template -void DoFRenumbering::Cuthill_McKee -(MGDoFHandler&, - const unsigned int, - const bool, - const std::vector&); - -template -void DoFRenumbering::cell_wise_dg -(MGDoFHandler&, - const unsigned int, - const std::vector::cell_iterator>&); -template -void DoFRenumbering::compute_cell_wise_dg -(std::vector&, std::vector&, - const MGDoFHandler&, const unsigned int, - const std::vector::cell_iterator>&); -template -void DoFRenumbering::compute_downstream_dg -(std::vector&, std::vector&, - const MGDoFHandler&, const unsigned int, - const Point&); + template + void random > + (DoFHandler &); + + template + void + compute_random > + (std::vector&, + const DoFHandler &); + + template + void sort_selected_dofs_back > + (DoFHandler &, + const std::vector &); + + template + void + compute_sort_selected_dofs_back > + (std::vector&, + const DoFHandler &, + const std::vector &); + + template + void random > + (hp::DoFHandler &); + + template + void + compute_random > + (std::vector&, + const hp::DoFHandler &); + + template + void sort_selected_dofs_back > + (hp::DoFHandler &, + const std::vector &); + + template + void + compute_sort_selected_dofs_back > + (std::vector&, + const hp::DoFHandler &, + const std::vector &); + + template + void subdomain_wise + (DoFHandler &); + + template + void + compute_subdomain_wise + (std::vector&, + const DoFHandler &); + + template + void Cuthill_McKee + (MGDoFHandler&, + const unsigned int, + const bool, + const std::vector&); + + template + void cell_wise_dg + (MGDoFHandler&, + const unsigned int, + const std::vector::cell_iterator>&); + template + void compute_cell_wise_dg + (std::vector&, std::vector&, + const MGDoFHandler&, const unsigned int, + const std::vector::cell_iterator>&); + template + void compute_downstream_dg + (std::vector&, std::vector&, + const MGDoFHandler&, const unsigned int, + const Point&); + +} // namespace DoFRenumbering DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index bed76caf08..a18014e31a 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -283,6 +283,21 @@ constraints individually.

deal.II

    +

  1. New: There is now a namespace DoFRenumbering::boost that contains the + implementation of interfaces to three reordering strategies provided by the + Boost Graph Library, namely DoFRenumbering::boost::Cuthill_McKee, + DoFRenumbering::boost::king_ordering, and DoFRenumbering::boost::minimum_degree. +
    + (WB 2008/03/07) +

    2. + +

  2. Changed: The DoFRenumbering class has been converted into a namespace. + Since all its members were static, this hasn't changed anything fundamental + and use of the functions therein should work in exactly the same way as before. +
    + (WB 2008/03/07) +

    3. +

  3. New: The new TriaObjectAccessor::set_all_boundary_indicators function does not only set the boundary indicator of the current object but of all that bound it as well. For example, in 3d, if TriaObjectAccessor::set_boundary_indicator @@ -293,8 +308,7 @@ constraints individually. (WB 2008/02/26)

    4. - -

  4. Improved: The tria pointer and the debug_output_grid function in GridIn +

  5. Improved: The tria pointer and the GridIn::debug_output_grid function are changed from private to protected to be accessible in derived classes.
    (RH 2008/02/24) -- 2.39.5