From: kanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d> Date: Mon, 30 Apr 2007 03:13:31 +0000 (+0000) Subject: add cellwise numbering for continuous methods X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=449f320a467cf2895a38c21a3df8cc50335c3f31;p=dealii-svn.git add cellwise numbering for continuous methods git-svn-id: https://svn.dealii.org/trunk@14646 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/dof_renumbering.h b/deal.II/deal.II/include/dofs/dof_renumbering.h index 00f9c250e0..eb7c547e06 100644 --- a/deal.II/deal.II/include/dofs/dof_renumbering.h +++ b/deal.II/deal.II/include/dofs/dof_renumbering.h @@ -171,20 +171,25 @@ DEAL_II_NAMESPACE_OPEN * block. * * - * <h3>Cell-wise numbering for Discontinuous Galerkin FEM</h3> + * <h3>Cell-wise numbering</h3> * - * 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. + * Given an ordered vector of cells, the function cell_wise() + * sorts the degrees of freedom such that degrees on earlier cells of + * this vector will occur before degrees on later cells. * - * 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. + * This rule produces a well-defined ordering for discontinuous Galerkin + * methods (FE_DGP, FE_DGQ). For continuous methods, we use the + * additional rule that each degree of freedom is ordered according to + * the first cell in the ordered vector it belongs to. * - * Given an ordered vector of cells, the function cell_wise_dg() - * accomplishes this. Inside the cells, the previous ordering will be - * preserved, so it may be useful to apply component_wise() first. + * Applications of this scheme are downstream() and + * clock_wise_dg(). The first orders the cells according to a + * downstream direction and then applies cell_wise(). + * + * @note For DG elements, the internal numbering in each cell remains + * unaffected. This cannot be guaranteed for continuous elements + * anymore, since degrees of freedom shared with an earlier cell will + * be accounted for by the other cell. * * * <h3>Random renumbering</h3> @@ -203,7 +208,7 @@ DEAL_II_NAMESPACE_OPEN * information on this. * * @ingroup dofs - * @author Wolfgang Bangerth, Guido Kanschat, 1998, 1999, 2000, 2004 + * @author Wolfgang Bangerth, Guido Kanschat, 1998, 1999, 2000, 2004, 2007 */ class DoFRenumbering { @@ -403,6 +408,16 @@ class DoFRenumbering */ template <class DH> static void + cell_wise (DH &dof_handler, + const std::vector<typename DH::cell_iterator> &cell_order); + + /** + * @deprecated Use cell_wise() instead. + * + * Cell-wise numbering only for DG. + */ + template <class DH> + static void cell_wise_dg (DH &dof_handler, const std::vector<typename DH::cell_iterator> &cell_order); @@ -423,6 +438,35 @@ class DoFRenumbering const std::vector<typename DH::cell_iterator> &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 <class DH> + static void + compute_cell_wise (std::vector<unsigned int>& renumbering, + std::vector<unsigned int>& inverse_renumbering, + const DH &dof_handler, + const std::vector<typename DH::cell_iterator> &cell_order); + + /** + * Cell-wise renumbering on one + * level. See the other function + * with the same name. + */ + template <int dim> + static void + cell_wise (MGDoFHandler<dim> &dof_handler, + const unsigned int level, + const std::vector<typename MGDoFHandler<dim>::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 @@ -452,6 +496,24 @@ class DoFRenumbering const std::vector<typename MGDoFHandler<dim>::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 <int dim> + static void + compute_cell_wise (std::vector<unsigned int>& renumbering, + std::vector<unsigned int>& inverse_renumbering, + const MGDoFHandler<dim>& dof_handler, + const unsigned int level, + const std::vector<typename MGDoFHandler<dim>::cell_iterator>& cell_order); + + /** * Cell-wise downstream numbering * with respect to a constant @@ -481,6 +543,14 @@ class DoFRenumbering */ template <class DH, int dim> static void + downstream (DH& dof_handler, + const Point<dim>& direction); + + /** + * @deprecated Use downstream() instead. + */ + template <class DH, int dim> + static void downstream_dg (DH& dof_handler, const Point<dim>& direction); @@ -493,6 +563,15 @@ class DoFRenumbering */ template <int dim> static void + downstream (MGDoFHandler<dim> &dof_handler, + const unsigned int level, + const Point<dim> &direction); + /** + * @deprecated Use downstream() + * instead. + */ + template <int dim> + static void downstream_dg (MGDoFHandler<dim> &dof_handler, const unsigned int level, const Point<dim> &direction); @@ -530,6 +609,17 @@ class DoFRenumbering */ template <class DH, int dim> static void + compute_downstream (std::vector<unsigned int>& new_dof_indices, + std::vector<unsigned int>& reverse, + const DH& dof_handler, + const Point<dim>& direction); + + /** + * @deprecated Use + * compute_downstream() instead + */ + template <class DH, int dim> + static void compute_downstream_dg (std::vector<unsigned int>& new_dof_indices, std::vector<unsigned int>& reverse, const DH& dof_handler, @@ -546,6 +636,18 @@ class DoFRenumbering */ template <int dim> static void + compute_downstream (std::vector<unsigned int>& new_dof_indices, + std::vector<unsigned int>& reverse, + const MGDoFHandler<dim>& dof_handler, + const unsigned int level, + const Point<dim>& direction); + + /** + * @deprecated Use + * compute_downstream() instead + */ + template <int dim> + static void compute_downstream_dg (std::vector<unsigned int>& new_dof_indices, std::vector<unsigned int>& reverse, const MGDoFHandler<dim>& dof_handler, diff --git a/deal.II/deal.II/source/dofs/dof_renumbering.cc b/deal.II/deal.II/source/dofs/dof_renumbering.cc index be8047068c..3a9f5ec507 100644 --- a/deal.II/deal.II/source/dofs/dof_renumbering.cc +++ b/deal.II/deal.II/source/dofs/dof_renumbering.cc @@ -889,7 +889,21 @@ DoFRenumbering::cell_wise_dg ( } +template <class DH> +void +DoFRenumbering::cell_wise ( + DH& dof, + const std::vector<typename DH::cell_iterator>& cells) +{ + std::vector<unsigned int> renumbering(dof.n_dofs()); + std::vector<unsigned int> reverse(dof.n_dofs()); + compute_cell_wise(renumbering, reverse, dof, cells); + + dof.renumber_dofs(renumbering); +} + +//TODO: Discuss if we cannot replace this function by the next template <class DH> void DoFRenumbering::compute_cell_wise_dg ( @@ -948,6 +962,71 @@ DoFRenumbering::compute_cell_wise_dg ( } +template <class DH> +void +DoFRenumbering::compute_cell_wise ( + std::vector<unsigned int>& new_indices, + std::vector<unsigned int>& reverse, + const DH& dof, + const typename std::vector<typename DH::cell_iterator>& 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(); + + // 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)); + + // For continuous elements, we must + // make sure, that each dof is + // reordered only once. + std::vector<bool> already_sorted(n_global_dofs, false); + std::vector<unsigned int> cell_dofs; + + unsigned int global_index = 0; + + typename std::vector<typename DH::cell_iterator>::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); + + (*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;i<n_cell_dofs;++i) + { + if (!already_sorted[cell_dofs[i]]) + { + already_sorted[cell_dofs[i]] = true; + reverse[global_index++] = cell_dofs[i]; + } + } + } + Assert(global_index == n_global_dofs, ExcRenumberingIncomplete()); + + for (unsigned int i=0;i<reverse.size(); ++i) + new_indices[reverse[i]] = i; +} + + template <int dim> void DoFRenumbering::cell_wise_dg ( MGDoFHandler<dim>& dof, @@ -958,7 +1037,21 @@ void DoFRenumbering::cell_wise_dg ( std::vector<unsigned int> reverse(dof.n_dofs(level)); compute_cell_wise_dg(renumbering, reverse, dof, level, cells); - dof.renumber_dofs(level, reverse); + dof.renumber_dofs(level, renumbering); +} + + +template <int dim> +void DoFRenumbering::cell_wise ( + MGDoFHandler<dim>& dof, + const unsigned int level, + const typename std::vector<typename MGDoFHandler<dim>::cell_iterator>& cells) +{ + std::vector<unsigned int> renumbering(dof.n_dofs(level)); + std::vector<unsigned int> reverse(dof.n_dofs(level)); + + compute_cell_wise(renumbering, reverse, dof, level, cells); + dof.renumber_dofs(level, renumbering); } @@ -1009,13 +1102,63 @@ void DoFRenumbering::compute_cell_wise_dg ( for (unsigned int i=0;i<n_cell_dofs;++i) { - new_order[global_index++] = cell_dofs[i]; + reverse[global_index++] = cell_dofs[i]; + } + } + Assert(global_index == n_global_dofs, ExcRenumberingIncomplete()); + + for (unsigned int i=0;i<new_order.size(); ++i) + new_order[reverse[i]] = i; +} + + + +template <int dim> +void DoFRenumbering::compute_cell_wise ( + std::vector<unsigned int>& new_order, + std::vector<unsigned int>& reverse, + const MGDoFHandler<dim>& dof, + const unsigned int level, + const typename std::vector<typename MGDoFHandler<dim>::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(); + + std::vector<bool> already_sorted(n_global_dofs, false); + std::vector<unsigned int> cell_dofs(n_cell_dofs); + + unsigned int global_index = 0; + + typename std::vector<typename MGDoFHandler<dim>::cell_iterator>::const_iterator cell; + + 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()); + + for (unsigned int i=0;i<n_cell_dofs;++i) + { + if (!already_sorted[cell_dofs[i]]) + { + already_sorted[cell_dofs[i]] = true; + reverse[global_index++] = cell_dofs[i]; + } } } Assert(global_index == n_global_dofs, ExcRenumberingIncomplete()); for (unsigned int i=0;i<new_order.size(); ++i) - reverse[new_order[i]] = i; + new_order[reverse[i]] = i; } @@ -1032,6 +1175,19 @@ DoFRenumbering::downstream_dg (DH& dof, const Point<dim>& direction) +template <class DH, int dim> +void +DoFRenumbering::downstream (DH& dof, const Point<dim>& direction) +{ + std::vector<unsigned int> renumbering(dof.n_dofs()); + std::vector<unsigned int> reverse(dof.n_dofs()); + compute_downstream(renumbering, reverse, dof, direction); + + dof.renumber_dofs(renumbering); +} + + + template <class DH, int dim> void DoFRenumbering::compute_downstream_dg ( @@ -1076,10 +1232,62 @@ DoFRenumbering::compute_downstream_dg ( } +template <class DH, int dim> +void +DoFRenumbering::compute_downstream ( + std::vector<unsigned int>& new_indices, + std::vector<unsigned int>& reverse, + const DH& dof, + const Point<dim>& direction) +{ + std::vector<typename DH::cell_iterator> + ordered_cells(dof.get_tria().n_active_cells()); + const CompareDownstream<dim> comparator(direction); + + typename DH::active_cell_iterator begin = dof.begin_active(); + typename DH::active_cell_iterator end = dof.end(); + + copy (begin, end, ordered_cells.begin()); + sort (ordered_cells.begin(), ordered_cells.end(), comparator); + + compute_cell_wise(new_indices, reverse, dof, ordered_cells); +} + + template <int dim> void DoFRenumbering::downstream_dg (MGDoFHandler<dim>& dof, const unsigned int level, const Point<dim>& direction) +{ + std::vector<unsigned int> renumbering(dof.n_dofs(level)); + std::vector<unsigned int> reverse(dof.n_dofs(level)); + compute_downstream_dg(renumbering, reverse, dof, level, direction); + + dof.renumber_dofs(level, renumbering); +} + + +template <int dim> +void DoFRenumbering::downstream (MGDoFHandler<dim>& dof, + const unsigned int level, + const Point<dim>& direction) +{ + std::vector<unsigned int> renumbering(dof.n_dofs(level)); + std::vector<unsigned int> reverse(dof.n_dofs(level)); + compute_downstream(renumbering, reverse, dof, level, direction); + + dof.renumber_dofs(level, renumbering); +} + + +template <int dim> +void +DoFRenumbering::compute_downstream_dg ( + std::vector<unsigned int>& new_indices, + std::vector<unsigned int>& reverse, + const MGDoFHandler<dim>& dof, + const unsigned int level, + const Point<dim>& direction) { std::vector<typename MGDoFHandler<dim>::cell_iterator> ordered_cells(dof.get_tria().n_cells(level)); @@ -1088,16 +1296,16 @@ void DoFRenumbering::downstream_dg (MGDoFHandler<dim>& dof, typename MGDoFHandler<dim>::cell_iterator begin = dof.begin(level); typename MGDoFHandler<dim>::cell_iterator end = dof.end(level); - std::copy (begin, end, ordered_cells.begin()); - std::sort (ordered_cells.begin(), ordered_cells.end(), comparator); - - cell_wise_dg (dof, level, ordered_cells); + copy (begin, end, ordered_cells.begin()); + sort (ordered_cells.begin(), ordered_cells.end(), comparator); + + compute_cell_wise_dg(new_indices, reverse, dof, level, ordered_cells); } template <int dim> void -DoFRenumbering::compute_downstream_dg ( +DoFRenumbering::compute_downstream ( std::vector<unsigned int>& new_indices, std::vector<unsigned int>& reverse, const MGDoFHandler<dim>& dof, @@ -1114,7 +1322,7 @@ DoFRenumbering::compute_downstream_dg ( copy (begin, end, ordered_cells.begin()); sort (ordered_cells.begin(), ordered_cells.end(), comparator); - compute_cell_wise_dg(new_indices, reverse, dof, level, ordered_cells); + compute_cell_wise(new_indices, reverse, dof, level, ordered_cells); } @@ -1403,43 +1611,67 @@ DoFRenumbering::compute_downstream_dg<DoFHandler<deal_II_dimension> > const DoFHandler<deal_II_dimension>&, const Point<deal_II_dimension>&); +template void +DoFRenumbering::cell_wise<DoFHandler<deal_II_dimension> > +(DoFHandler<deal_II_dimension>&, + const std::vector<DoFHandler<deal_II_dimension>::cell_iterator>&); + +template void +DoFRenumbering::compute_cell_wise<DoFHandler<deal_II_dimension> > +(std::vector<unsigned int>&, std::vector<unsigned int>&, + const DoFHandler<deal_II_dimension>&, + const std::vector<DoFHandler<deal_II_dimension>::cell_iterator>&); + +template void +DoFRenumbering::downstream<DoFHandler<deal_II_dimension> > +(DoFHandler<deal_II_dimension>&, const Point<deal_II_dimension>&); + +template void +DoFRenumbering::compute_downstream<DoFHandler<deal_II_dimension> > +(std::vector<unsigned int>&,std::vector<unsigned int>&, + const DoFHandler<deal_II_dimension>&, const Point<deal_II_dimension>&); + template void DoFRenumbering::clockwise_dg<DoFHandler<deal_II_dimension> > -(DoFHandler<deal_II_dimension>&, - const Point<deal_II_dimension>&, bool); +(DoFHandler<deal_II_dimension>&, const Point<deal_II_dimension>&, bool); template void DoFRenumbering::compute_clockwise_dg<DoFHandler<deal_II_dimension> > -(std::vector<unsigned int>&, - const DoFHandler<deal_II_dimension>&, - const Point<deal_II_dimension>&, - const bool); +(std::vector<unsigned int>&, const DoFHandler<deal_II_dimension>&, + const Point<deal_II_dimension>&, const bool); -// DG renumbering for hp::DoFHandler +// Renumbering for hp::DoFHandler -template -void +template void DoFRenumbering::cell_wise_dg<hp::DoFHandler<deal_II_dimension> > (hp::DoFHandler<deal_II_dimension>&, const std::vector<hp::DoFHandler<deal_II_dimension>::cell_iterator>&); -template -void +template void DoFRenumbering::compute_cell_wise_dg<hp::DoFHandler<deal_II_dimension> > (std::vector<unsigned int>&, std::vector<unsigned int>&, const hp::DoFHandler<deal_II_dimension>&, const std::vector<hp::DoFHandler<deal_II_dimension>::cell_iterator>&); -template -void +template void +DoFRenumbering::cell_wise<hp::DoFHandler<deal_II_dimension> > +(hp::DoFHandler<deal_II_dimension>&, + const std::vector<hp::DoFHandler<deal_II_dimension>::cell_iterator>&); + +template void +DoFRenumbering::compute_cell_wise<hp::DoFHandler<deal_II_dimension> > +(std::vector<unsigned int>&, std::vector<unsigned int>&, + const hp::DoFHandler<deal_II_dimension>&, + const std::vector<hp::DoFHandler<deal_II_dimension>::cell_iterator>&); + +template void DoFRenumbering::downstream_dg<hp::DoFHandler<deal_II_dimension> > (hp::DoFHandler<deal_II_dimension>&, const Point<deal_II_dimension>&); -template -void +template void DoFRenumbering::compute_downstream_dg<hp::DoFHandler<deal_II_dimension> > (std::vector<unsigned int>&,std::vector<unsigned int>&, const hp::DoFHandler<deal_II_dimension>&, @@ -1451,6 +1683,17 @@ DoFRenumbering::compute_downstream_dg<hp::DoFHandler<deal_II_dimension> > const hp::DoFHandler<deal_II_dimension>&, const Point<deal_II_dimension>&); +template void +DoFRenumbering::downstream<hp::DoFHandler<deal_II_dimension> > +(hp::DoFHandler<deal_II_dimension>&, + const Point<deal_II_dimension>&); + +template void +DoFRenumbering::compute_downstream<hp::DoFHandler<deal_II_dimension> > +(std::vector<unsigned int>&,std::vector<unsigned int>&, + const hp::DoFHandler<deal_II_dimension>&, + const Point<deal_II_dimension>&); + template void DoFRenumbering::clockwise_dg<hp::DoFHandler<deal_II_dimension> > @@ -1469,8 +1712,12 @@ DoFRenumbering::compute_clockwise_dg<hp::DoFHandler<deal_II_dimension> > template void DoFRenumbering::downstream_dg<deal_II_dimension> -(MGDoFHandler<deal_II_dimension>&, - const unsigned int, +(MGDoFHandler<deal_II_dimension>&, const unsigned int, + const Point<deal_II_dimension>&); + +template +void DoFRenumbering::downstream<deal_II_dimension> +(MGDoFHandler<deal_II_dimension>&, const unsigned int, const Point<deal_II_dimension>&); template