const bool use_constraints = false,
const std::vector<unsigned int> &starting_indices = std::vector<unsigned int>());
+ /**
+ * Computes the renumbering
+ * vector needed by the
+ * @p{Cuthill_McKee}
+ * function. Does not perform the
+ * renumbering on the
+ * @p{DoFHandler} dofs but
+ * returns the renumbering
+ * vector.
+ */
+ template <int dim>
+ static std::vector<unsigned int>
+ compute_Cuthill_McKee (DoFHandler<dim> &dof_handler,
+ const bool reversed_numbering = false,
+ const bool use_constraints = false,
+ const std::vector<unsigned int> &starting_indices = std::vector<unsigned int>());
+
/**
* Renumber the degrees of freedom
* according to the Cuthill-McKee method,
*/
template <int dim>
static void
- component_wise (DoFHandler<dim> &dof_handler,
+ component_wise (DoFHandler<dim> &dof_handler,
const std::vector<unsigned int> &component_order = std::vector<unsigned int>());
+ /**
+ * Computes the renumbering
+ * vector needed by the
+ * @p{component_wise}
+ * function. Does not perform the
+ * renumbering on the
+ * @p{DoFHandler} dofs but
+ * returns the renumbering
+ * vector.
+ */
+ template <int dim>
+ static std::vector<unsigned int>
+ compute_component_wise (DoFHandler<dim> &dof_handler,
+ const std::vector<unsigned int> &component_order = std::vector<unsigned int>());
+
/**
* Sort the degrees of freedom by
* component. It does the same
*/
template <int dim>
static void
- cell_wise_dg (DoFHandler<dim> &dof_handler,
+ cell_wise_dg (DoFHandler<dim> &dof_handler,
const std::vector<typename DoFHandler<dim>::cell_iterator> &cell_order);
-
+
+ /**
+ * Computes the renumbering
+ * vector needed by the
+ * @p{cell_wise_dg}
+ * function. Does not perform the
+ * renumbering on the
+ * @p{DoFHandler} dofs but
+ * returns the renumbering
+ * vector.
+ */
+ template <int dim>
+ static std::vector<unsigned int>
+ compute_cell_wise_dg (DoFHandler<dim> &dof_handler,
+ const std::vector<typename DoFHandler<dim>::cell_iterator> &cell_order);
/**
* Cell-wise renumbering on one
downstream_dg (DoFHandler<dim> &dof_handler,
const Point<dim> &direction);
+ /**
+ * Computes the renumbering
+ * vector needed by the
+ * @p{downstream_dg}
+ * function. Does not perform the
+ * renumbering on the
+ * @p{DoFHandler} dofs but
+ * returns the renumbering
+ * vector.
+ */
+ template <int dim>
+ static std::vector<unsigned int>
+ compute_downstream_dg (DoFHandler<dim> &dof_handler,
+ const Point<dim> &direction);
+
/**
* Cell-wise downstream numbering
* with respect to a constant
sort_selected_dofs_back (DoFHandler<dim> &dof_handler,
const std::vector<bool> &selected_dofs);
+ /**
+ * Computes the renumbering
+ * vector needed by the
+ * @p{sort_selected_dofs_back}
+ * function. Does not perform the
+ * renumbering on the
+ * @p{DoFHandler} dofs but
+ * returns the renumbering
+ * vector.
+ */
+ template <int dim>
+ static std::vector<unsigned int>
+ compute_sort_selected_dofs_back (DoFHandler<dim> &dof_handler,
+ const std::vector<bool> &selected_dofs);
+
/**
* Renumber the degrees of
* freedom in a random way.
template <int dim>
static void
random (DoFHandler<dim> &dof_handler);
-
+
+ /**
+ * Computes the renumbering
+ * vector needed by the
+ * @p{random} function. Does not
+ * perform the renumbering on the
+ * @p{DoFHandler} dofs but
+ * returns the renumbering
+ * vector.
+ */
+ template <int dim>
+ static std::vector<unsigned int>
+ compute_random (DoFHandler<dim> &dof_handler);
+
/**
* Exception
*/
#endif
+
template <int dim>
-void DoFRenumbering::Cuthill_McKee (DoFHandler<dim> &dof_handler,
- const bool reversed_numbering,
- const bool use_constraints,
- const std::vector<unsigned int> &starting_indices)
+void
+DoFRenumbering::Cuthill_McKee (DoFHandler<dim> &dof_handler,
+ const bool reversed_numbering,
+ const bool use_constraints,
+ const std::vector<unsigned int> &starting_indices)
+{
+ std::vector<unsigned int> renumbering
+ =compute_Cuthill_McKee(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 <int dim>
+std::vector<unsigned int>
+DoFRenumbering::compute_Cuthill_McKee (DoFHandler<dim> &dof_handler,
+ const bool reversed_numbering,
+ const bool use_constraints,
+ const std::vector<unsigned int> &starting_indices)
{
// make the connection graph
SparsityPattern sparsity (dof_handler.n_dofs(),
i!=new_number.end(); ++i)
*i = n_dofs-*i-1;
- // actually perform renumbering;
- // this is dimension specific and
- // thus needs an own function
- dof_handler.renumber_dofs (new_number);
+ return new_number;
}
void
DoFRenumbering::component_wise (DoFHandler<dim> &dof_handler,
const std::vector<unsigned int> &component_order_arg)
+{
+ std::vector<unsigned int> renumbering
+ =compute_component_wise(dof_handler, component_order_arg);
+
+ if (renumbering.size()!=0)
+ dof_handler.renumber_dofs (renumbering);
+}
+
+
+
+template <int dim>
+std::vector<unsigned int>
+DoFRenumbering::compute_component_wise (DoFHandler<dim> &dof_handler,
+ const std::vector<unsigned int> &component_order_arg)
{
const unsigned int dofs_per_cell = dof_handler.get_fe().dofs_per_cell;
const FiniteElement<dim> &fe = dof_handler.get_fe();
+ std::vector<unsigned int> new_indices (dof_handler.n_dofs(),
+ DoFHandler<dim>::invalid_dof_index);
+
// do nothing if the FE has only
// one component
if (fe.n_components() == 1)
- return;
+ {
+ new_indices.resize(0);
+ return new_indices;
+ }
std::vector<unsigned int> component_order (component_order_arg);
if (component_order.size() == 0)
// components in the order the user
// desired to see
unsigned int next_free_index = 0;
- std::vector<unsigned int> new_indices (dof_handler.n_dofs(),
- DoFHandler<dim>::invalid_dof_index);
for (unsigned int c=0; c<fe.n_components(); ++c)
{
const unsigned int component = component_order[c];
Assert (next_free_index == dof_handler.n_dofs(),
ExcInternalError());
- dof_handler.renumber_dofs (new_indices);
+ return new_indices;
}
void
DoFRenumbering::sort_selected_dofs_back (DoFHandler<dim> &dof_handler,
const std::vector<bool> &selected_dofs)
+{
+ std::vector<unsigned int> renumbering
+ =compute_sort_selected_dofs_back(dof_handler, selected_dofs);
+
+ dof_handler.renumber_dofs(renumbering);
+}
+
+
+
+template <int dim>
+std::vector<unsigned int>
+DoFRenumbering::compute_sort_selected_dofs_back (DoFHandler<dim> &dof_handler,
+ const std::vector<bool> &selected_dofs)
{
const unsigned int n_dofs = dof_handler.n_dofs();
Assert (selected_dofs.size() == n_dofs,
Assert (next_unselected == n_selected_dofs, ExcInternalError());
Assert (next_selected == n_dofs, ExcInternalError());
- // now perform the renumbering
- dof_handler.renumber_dofs (new_dof_indices);
+ return new_dof_indices;
}
+
template <int dim>
void
-DoFRenumbering::cell_wise_dg (DoFHandler<dim>& dof,
- const typename std::vector<typename DoFHandler<dim>::cell_iterator>& cells)
+DoFRenumbering::cell_wise_dg (
+ DoFHandler<dim>& dof,
+ const typename std::vector<typename DoFHandler<dim>::cell_iterator>& cells)
+{
+ std::vector<unsigned int> renumbering
+ =compute_cell_wise_dg(dof, cells);
+
+ dof.renumber_dofs(renumbering);
+}
+
+
+
+template <int dim>
+std::vector<unsigned int>
+DoFRenumbering::compute_cell_wise_dg (
+ DoFHandler<dim>& dof,
+ const typename std::vector<typename DoFHandler<dim>::cell_iterator>& cells)
{
Assert(cells.size() == dof.get_tria().n_active_cells(),
ExcDimensionMismatch(cells.size(),
for (unsigned int i=0;i<new_order.size(); ++i)
reverse[new_order[i]] = i;
- dof.renumber_dofs(reverse);
+ return reverse;
}
#ifdef ENABLE_MULTIGRID
template <int dim>
-void
-DoFRenumbering::cell_wise_dg (MGDoFHandler<dim>& dof,
- const unsigned int level,
- const typename std::vector<typename MGDoFHandler<dim>::cell_iterator>& cells)
+void DoFRenumbering::cell_wise_dg (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(),
}
};
-
+
template <int dim>
void
DoFRenumbering::downstream_dg (DoFHandler<dim>& dof,
const Point<dim>& direction)
+{
+ std::vector<unsigned int> renumbering
+ =compute_downstream_dg(dof, direction);
+
+ dof.renumber_dofs(renumbering);
+}
+
+
+
+template <int dim>
+std::vector<unsigned int>
+DoFRenumbering::compute_downstream_dg (DoFHandler<dim>& dof,
+ const Point<dim>& direction)
{
std::vector<typename DoFHandler<dim>::cell_iterator>
ordered_cells(dof.get_tria().n_active_cells());
copy (begin, end, ordered_cells.begin());
sort (ordered_cells.begin(), ordered_cells.end(), comparator);
- cell_wise_dg(dof, ordered_cells);
+ return compute_cell_wise_dg(dof, ordered_cells);
}
#ifdef ENABLE_MULTIGRID
template <int dim>
-void
-DoFRenumbering::downstream_dg (MGDoFHandler<dim>& dof,
- const unsigned int level,
- const Point<dim>& direction)
+void DoFRenumbering::downstream_dg (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));
#endif
+
template <int dim>
void
DoFRenumbering::random (DoFHandler<dim> &dof_handler)
+{
+ std::vector<unsigned int> renumbering
+ =compute_random(dof_handler);
+
+ dof_handler.renumber_dofs(renumbering);
+}
+
+
+template <int dim>
+std::vector<unsigned int>
+DoFRenumbering::compute_random (DoFHandler<dim> &dof_handler)
{
const unsigned int n_dofs = dof_handler.n_dofs();
new_indices[i] = i;
std::random_shuffle (new_indices.begin(), new_indices.end());
- dof_handler.renumber_dofs(new_indices);
+ return new_indices;
}
const bool,
const std::vector<unsigned int>&);
+template
+std::vector<unsigned int>
+DoFRenumbering::compute_Cuthill_McKee<deal_II_dimension>
+(DoFHandler<deal_II_dimension>&,
+ const bool,
+ const bool,
+ const std::vector<unsigned int>&);
+
template
void DoFRenumbering::component_wise<deal_II_dimension>
(DoFHandler<deal_II_dimension>&,
const std::vector<unsigned int>&);
+template
+std::vector<unsigned int>
+DoFRenumbering::compute_component_wise<deal_II_dimension>
+(DoFHandler<deal_II_dimension>&,
+ const std::vector<unsigned int>&);
+
template
void DoFRenumbering::component_wise<deal_II_dimension>
(MGDoFHandler<deal_II_dimension>&,
template
-void DoFRenumbering::cell_wise_dg<deal_II_dimension>
+void
+DoFRenumbering::cell_wise_dg<deal_II_dimension>
(DoFHandler<deal_II_dimension>&,
const std::vector<DoFHandler<deal_II_dimension>::cell_iterator>&);
template
-void DoFRenumbering::downstream_dg<deal_II_dimension>
+std::vector<unsigned int>
+DoFRenumbering::compute_cell_wise_dg<deal_II_dimension>
+(DoFHandler<deal_II_dimension>&,
+ const std::vector<DoFHandler<deal_II_dimension>::cell_iterator>&);
+
+template
+void
+DoFRenumbering::downstream_dg<deal_II_dimension>
+(DoFHandler<deal_II_dimension>&,
+ const Point<deal_II_dimension>&);
+
+template
+std::vector<unsigned int>
+DoFRenumbering::compute_downstream_dg<deal_II_dimension>
(DoFHandler<deal_II_dimension>&,
const Point<deal_II_dimension>&);
(DoFHandler<deal_II_dimension> &,
const std::vector<bool> &);
+template
+std::vector<unsigned int>
+DoFRenumbering::compute_sort_selected_dofs_back<deal_II_dimension>
+(DoFHandler<deal_II_dimension> &,
+ const std::vector<bool> &);
+
template
void DoFRenumbering::random<deal_II_dimension>
(DoFHandler<deal_II_dimension> &);
+template
+std::vector<unsigned int>
+DoFRenumbering::compute_random<deal_II_dimension>
+(DoFHandler<deal_II_dimension> &);
+
#ifdef ENABLE_MULTIGRID
template
void DoFRenumbering::Cuthill_McKee<deal_II_dimension>
const unsigned int,
const bool,
const std::vector<unsigned int>&);
+
template
void DoFRenumbering::cell_wise_dg<deal_II_dimension>
(MGDoFHandler<deal_II_dimension>&,