From cfa0590af3628567140d8ce20a977d2a9b605518 Mon Sep 17 00:00:00 2001 From: guido Date: Thu, 13 Jul 2000 14:33:43 +0000 Subject: [PATCH] cell_wise numbering in DoFRenumbering git-svn-id: https://svn.dealii.org/trunk@3157 0785d39b-7218-0410-832d-ea1e28bc413d --- .../include/numerics/dof_renumbering.h | 53 +++++++++++-- .../deal.II/source/dofs/dof_renumbering.cc | 75 +++++++++++++++---- 2 files changed, 108 insertions(+), 20 deletions(-) diff --git a/deal.II/deal.II/include/numerics/dof_renumbering.h b/deal.II/deal.II/include/numerics/dof_renumbering.h index 7b95b83d52..6d3d61e495 100644 --- a/deal.II/deal.II/include/numerics/dof_renumbering.h +++ b/deal.II/deal.II/include/numerics/dof_renumbering.h @@ -139,7 +139,7 @@ * may be difficult, however, and in many cases will not justify the effort. * * - * @sect2{Componentwise numbering} + * @sect2{Component-wise numbering} * * For finite elements composed of several base elements using the @p{FESystem} * class, or for elements which provide several components themselves, it @@ -154,6 +154,20 @@ * algorithm and afterwards renumbering component-wise. This will bring out the * matrix structure and additionally have a good numbering within each block. * + * @sect2{Cell-wise numbering for Discontinuous Galerkin FEM} + * + * One advantage of DGFEM is the fact, that it yields invertible + * blocks on the diagonal of the global matrix. these blocks are in + * fact the cell matrices and matrix elements outside these blocks are + * due to fluxes. + * + * Still, it may be necessary to apply a downstream numbering of the + * degrees of freedom. This renumbering may only exchange whole blocks + * and must not destroy the block structure. + * + * Given an ordered vector of cells, the function @p{cell_wise} + * accomplishes this. Inside the cells, the previous ordering will be + * preserved, so it may be useful to apply @component_wise} first. * * @sect2{Multigrid DoF numbering} * @@ -161,7 +175,7 @@ * to the actual function declarations to get more information on this. * * - * @author Wolfgang Bangerth, 1998, 1999 + * @author Wolfgang Bangerth, Guido Kanschat, 1998, 1999, 2000 */ class DoFRenumbering { @@ -245,6 +259,31 @@ class DoFRenumbering component_wise (DoFHandler &dof_handler, const vector &component_order = vector()); + /** + * Cell-wise renumbering for DG + * elements. This function takes + * the ordered set of cells in + * @p{cell_order}, and makes sure + * that all degrees of freedom in + * a cell with higher index are + * behind all degrees of freedom + * of a cell with lower + * index. The order inside a cell + * bloock will be the same as + * before this renumbering. + * + * 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< int dim> + static void + cell_wise (DoFHandler &dof_handler, + const vector > &cell_order); + + /** * Sort those degrees of freedom * which are tagged with @p{true} @@ -271,13 +310,13 @@ class DoFRenumbering * Exception */ DeclException0 (ExcInvalidComponentOrder); + /** - * Exception + * The function is only + * implemented for Discontinuous + * Galerkin Finite elements. */ - DeclException2 (ExcInvalidArraySize, - int, int, - << "The array has size " << arg1 << " but should have " - << arg2 << "."); + DeclException0 (ExcNotDGFEM); }; diff --git a/deal.II/deal.II/source/dofs/dof_renumbering.cc b/deal.II/deal.II/source/dofs/dof_renumbering.cc index a7ad38139c..ae82208b61 100644 --- a/deal.II/deal.II/source/dofs/dof_renumbering.cc +++ b/deal.II/deal.II/source/dofs/dof_renumbering.cc @@ -15,6 +15,7 @@ #include #include #include +#include #include #include #include @@ -505,7 +506,7 @@ DoFRenumbering::sort_selected_dofs_back (const vector &selected_dofs, { const unsigned int n_dofs = dof_handler.n_dofs(); Assert (selected_dofs.size() == n_dofs, - ExcInvalidArraySize (selected_dofs.size(), n_dofs)); + ExcDimensionMismatch (selected_dofs.size(), n_dofs)); // re-sort the dofs according to // their selection state @@ -535,27 +536,75 @@ DoFRenumbering::sort_selected_dofs_back (const vector &selected_dofs, }; +template +void +DoFRenumbering::cell_wise (DoFHandler& dof, + const vector >& cells) +{ + Assert(cells.size() == dof.get_tria().n_active_cells(), + ExcDimensionMismatch(cells.size(), + dof.get_tria().n_active_cells())); + 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()); + } + + unsigned int n_global_dofs = dof.n_dofs(); + unsigned int n_cell_dofs = dof.get_fe().n_dofs_per_cell(); + + vector new_order(n_global_dofs); + vector cell_dofs(n_cell_dofs); + + unsigned int global_index = 0; + + typename vector >::const_iterator cell; + + for(cell = cells.begin(); cell != cells.end(); ++cell) + { + cell->get_dof_indices(cell_dofs); + sort(cell_dofs.begin(), cell_dofs.end()); + + for (unsigned int i=0;i &dof_handler, - const bool reversed_numbering, - const bool use_constraints, - const vector &starting_indices); +void DoFRenumbering::Cuthill_McKee (DoFHandler&, + const bool, + const bool, + const vector&); + +template +void DoFRenumbering::Cuthill_McKee (MGDoFHandler&, + const unsigned int, + const bool, + const vector&); template -void DoFRenumbering::Cuthill_McKee (MGDoFHandler &dof_handler, - const unsigned int level, - const bool reversed_numbering, - const vector &starting_indices); +void DoFRenumbering::component_wise (DoFHandler&, + const vector&); + template -void DoFRenumbering::component_wise (DoFHandler &dof_handler, - const vector &component_order_arg); +void DoFRenumbering::cell_wise (DoFHandler&, + const vector >&); + template -void DoFRenumbering::sort_selected_dofs_back (const vector &selected_dofs, - DoFHandler &dof_handler); +void DoFRenumbering::sort_selected_dofs_back (const vector &, + DoFHandler &); -- 2.39.5