// $Id$
// Version: $Name$
//
-// Copyright (C) 2003, 2004, 2005, 2006 by the deal.II authors
+// Copyright (C) 2003, 2004, 2005, 2006, 2007 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
*/
template <class DH>
static void
- compute_cell_wise_dg (std::vector<unsigned int> &new_dof_indices,
+ compute_cell_wise_dg (std::vector<unsigned int>& renumbering,
+ std::vector<unsigned int>& inverse_renumbering,
const DH &dof_handler,
const std::vector<typename DH::cell_iterator> &cell_order);
const unsigned int level,
const std::vector<typename MGDoFHandler<dim>::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 <int dim>
+ static void
+ compute_cell_wise_dg (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
const Point<dim> &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
const DH& dof_handler,
const Point<dim>& 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 <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,
+ const Point<dim>& 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 <int dim>
+ static void
+ compute_downstream_dg (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);
+
//TODO:[GK] Documentation!
/**
* Cell-wise clockwise numbering.
const std::vector<typename DH::cell_iterator>& cells)
{
std::vector<unsigned int> renumbering(dof.n_dofs());
- compute_cell_wise_dg(renumbering, dof, cells);
+ std::vector<unsigned int> reverse(dof.n_dofs());
+ compute_cell_wise_dg(renumbering, reverse, dof, cells);
dof.renumber_dofs(renumbering);
}
void
DoFRenumbering::compute_cell_wise_dg (
std::vector<unsigned int>& new_indices,
+ std::vector<unsigned int>& reverse,
const DH& dof,
const typename std::vector<typename DH::cell_iterator>& cells)
{
Assert(new_indices.size() == n_global_dofs,
ExcDimensionMismatch(new_indices.size(), n_global_dofs));
- std::vector<unsigned int> reverse(new_indices.size());
-
+ Assert(reverse.size() == n_global_dofs,
+ ExcDimensionMismatch(reverse.size(), n_global_dofs));
+
std::vector<unsigned int> cell_dofs;
unsigned int global_index = 0;
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_dg(renumbering, reverse, dof, level, cells);
+ dof.renumber_dofs(level, reverse);
+}
+
+
+template <int dim>
+void DoFRenumbering::compute_cell_wise_dg (
+ 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(),
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)));
+
unsigned int n_global_dofs = dof.n_dofs(level);
unsigned int n_cell_dofs = dof.get_fe().n_dofs_per_cell();
-
- std::vector<unsigned int> new_order(n_global_dofs);
+
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;
}
Assert(global_index == n_global_dofs, ExcRenumberingIncomplete());
- std::vector<unsigned int> reverse(new_order.size());
-
for (unsigned int i=0;i<new_order.size(); ++i)
- reverse[new_order[i]] = i;
-
- dof.renumber_dofs(level, reverse);
+ reverse[new_order[i]] = i;
}
copy (begin, end, ordered_cells.begin());
sort (ordered_cells.begin(), ordered_cells.end(), comparator);
- compute_cell_wise_dg(new_indices, dof, ordered_cells);
+ std::vector<unsigned int> reverse(new_indices.size());
+ compute_cell_wise_dg(new_indices, reverse, dof, ordered_cells);
+}
+
+
+template <class DH, int dim>
+void
+DoFRenumbering::compute_downstream_dg (
+ 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_dg(new_indices, reverse, dof, ordered_cells);
}
}
+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));
+ const CompareDownstream<dim> comparator(direction);
+
+ typename MGDoFHandler<dim>::cell_iterator begin = dof.begin(level);
+ typename MGDoFHandler<dim>::cell_iterator end = dof.end(level);
+
+ 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);
+}
+
+
/**
* Provide comparator for DoFCellAccessors
copy (begin, end, ordered_cells.begin());
sort (ordered_cells.begin(), ordered_cells.end(), comparator);
- compute_cell_wise_dg(new_indices, dof, ordered_cells);
+ std::vector<unsigned int> reverse(new_indices.size());
+ compute_cell_wise_dg(new_indices, reverse, dof, ordered_cells);
}
template
void
DoFRenumbering::compute_cell_wise_dg<DoFHandler<deal_II_dimension> >
-(std::vector<unsigned int>&,
+(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::compute_downstream_dg<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::compute_downstream_dg<DoFHandler<deal_II_dimension> >
(std::vector<unsigned int>&,
const DoFHandler<deal_II_dimension>&,
const Point<deal_II_dimension>&);
template
void
DoFRenumbering::compute_cell_wise_dg<hp::DoFHandler<deal_II_dimension> >
-(std::vector<unsigned int>&,
+(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::compute_downstream_dg<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::compute_downstream_dg<hp::DoFHandler<deal_II_dimension> >
(std::vector<unsigned int>&,
const hp::DoFHandler<deal_II_dimension>&,
const Point<deal_II_dimension>&);
(MGDoFHandler<deal_II_dimension>&,
const unsigned int,
const std::vector<MGDoFHandler<deal_II_dimension>::cell_iterator>&);
+template
+void DoFRenumbering::compute_cell_wise_dg<deal_II_dimension>
+(std::vector<unsigned int>&, std::vector<unsigned int>&,
+ const MGDoFHandler<deal_II_dimension>&, const unsigned int,
+ const std::vector<MGDoFHandler<deal_II_dimension>::cell_iterator>&);
+template
+void DoFRenumbering::compute_downstream_dg<deal_II_dimension>
+(std::vector<unsigned int>&, std::vector<unsigned int>&,
+ const MGDoFHandler<deal_II_dimension>&, const unsigned int,
+ const Point<deal_II_dimension>&);
DEAL_II_NAMESPACE_CLOSE