* 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>
* 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
{
*/
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);
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
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
*/
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);
*/
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);
*/
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,
*/
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,
}
+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 (
}
+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,
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);
}
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;
}
+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 (
}
+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));
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,
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);
}
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>&,
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> >
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