* Extract the indices of the
* degrees of freedom belonging
* to certain vector components
- * of a vector-valued finite
- * element. The bit vector
- * @p select defines, which
- * components of an
- * FESystem are to be
- * extracted from the
- * DoFHandler @p dof. The
- * entries in @p selected_dofs
- * corresponding to degrees of
- * freedom belonging to these
- * components are then flagged
- * @p true, while all others are
- * set to @p false.
+ * or blocks (if the last
+ * argument is <tt>true</tt>) of
+ * a vector-valued finite
+ * element. The bit vector @p
+ * select defines, which
+ * components or blocks of an
+ * FESystem are to be extracted
+ * from the DoFHandler @p
+ * dof. The entries in @p
+ * selected_dofs corresponding to
+ * degrees of freedom belonging
+ * to these components are then
+ * flagged @p true, while all
+ * others are set to @p false.
*
- * The size of
- * @p component_select shall
+ * The size of @p select must
* equal the number of components
- * in the finite element used by
- * @p dof. The size of
- * @p selected_dofs shall equal
- * <tt>dof_handler.n_dofs()</tt>. Previous
+ * or blocks in the FiniteElement
+ * used by @p dof, depending on
+ * the argument
+ * <tt>blocks</tt>. The size of
+ * @p selected_dofs must equal
+ * DoFHandler::n_dofs(). Previous
* contents of this array are
* overwritten.
*
* of its shape functions are
* non-zero in more than one
* vector component (which holds,
- * for example, for Nedelec or
- * Raviart-Thomas elements), then
+ * for example, for FE_Nedelec or
+ * FE_RaviartThomas elements), then
* shape functions cannot be
* associated with a single
* vector component. In this
template <int dim>
static void
extract_dofs (const DoFHandler<dim> &dof_handler,
- const std::vector<bool> &component_select,
- std::vector<bool> &selected_dofs);
+ const std::vector<bool> &select,
+ std::vector<bool> &selected_dofs,
+ bool blocks = false);
/**
* Do the same thing as
- * @p extract_dofs for one level
+ * extract_dofs() for one level
* of a multi-grid DoF numbering.
*/
template <int dim>
extract_level_dofs (const unsigned int level,
const MGDoFHandler<dim> &dof,
const std::vector<bool> &select,
- std::vector<bool> &selected_dofs);
+ std::vector<bool> &selected_dofs,
+ bool blocks = false);
/**
* Extract all degrees of freedom
template <int dim>
void
-DoFTools::extract_dofs (const DoFHandler<dim> &dof,
- const std::vector<bool> &component_select,
- std::vector<bool> &selected_dofs)
+DoFTools::extract_dofs (
+ const DoFHandler<dim> &dof,
+ const std::vector<bool> &component_select,
+ std::vector<bool> &selected_dofs,
+ bool count_by_blocks)
{
- Assert(component_select.size() == n_components(dof),
- ExcDimensionMismatch(component_select.size(), n_components(dof)));
+ const FiniteElement<dim> &fe = dof.get_fe();
+
+ if (count_by_blocks)
+ {
+ Assert(component_select.size() == fe.n_blocks(),
+ ExcDimensionMismatch(component_select.size(), fe.n_blocks()));
+ }
+ else
+ {
+ Assert(component_select.size() == n_components(dof),
+ ExcDimensionMismatch(component_select.size(), n_components(dof)));
+ }
+
Assert(selected_dofs.size() == dof.n_dofs(),
ExcDimensionMismatch(selected_dofs.size(), dof.n_dofs()));
};
- const FiniteElement<dim> &fe = dof.get_fe();
-
// preset all values by false
std::fill_n (selected_dofs.begin(), dof.n_dofs(), false);
// something interesting or not
std::vector<bool> local_selected_dofs (fe.dofs_per_cell, false);
for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
- if (fe.is_primitive(i))
+ if (count_by_blocks)
local_selected_dofs[i]
- = component_select[fe.system_to_component_index(i).first];
+ = component_select[fe.system_to_block_index(i).first];
else
- // if this shape function is
- // not primitive, then we have
- // to work harder. we have to
- // find out whether _any_ of
- // the vector components of
- // this element is selected or
- // not
- //
- // to do so, get the first and
- // last vector components of
- // the base element to which
- // the local dof with index i
- // belongs
- {
- unsigned int first_comp = 0;
- const unsigned int this_base = fe.system_to_base_index(i).first.first;
- const unsigned int this_multiplicity
- = fe.system_to_base_index(i).first.second;
-
- for (unsigned int b=0; b<this_base; ++b)
- first_comp += fe.base_element(b).n_components() *
- fe.element_multiplicity(b);
- for (unsigned int m=0; m<this_multiplicity; ++m)
- first_comp += fe.base_element(this_base).n_components();
- const unsigned int end_comp = first_comp +
- fe.base_element(this_base).n_components();
-
- Assert (first_comp < fe.n_components(), ExcInternalError());
- Assert (end_comp <= fe.n_components(), ExcInternalError());
-
- // now check whether any of
- // the components in between
- // is set
- for (unsigned int c=first_comp; c<end_comp; ++c)
- if (component_select[c] == true)
- {
- local_selected_dofs[i] = true;
- break;
- };
- };
+ if (fe.is_primitive(i))
+ local_selected_dofs[i]
+ = component_select[fe.system_to_component_index(i).first];
+ else
+ // if this shape function is
+ // not primitive, then we have
+ // to work harder. we have to
+ // find out whether _any_ of
+ // the vector components of
+ // this element is selected or
+ // not
+ //
+ // to do so, get the first and
+ // last vector components of
+ // the base element to which
+ // the local dof with index i
+ // belongs
+ {
+ unsigned int first_comp = 0;
+ const unsigned int this_base = fe.system_to_base_index(i).first.first;
+ const unsigned int this_multiplicity
+ = fe.system_to_base_index(i).first.second;
+
+ for (unsigned int b=0; b<this_base; ++b)
+ first_comp += fe.base_element(b).n_components() *
+ fe.element_multiplicity(b);
+ for (unsigned int m=0; m<this_multiplicity; ++m)
+ first_comp += fe.base_element(this_base).n_components();
+ const unsigned int end_comp = first_comp +
+ fe.base_element(this_base).n_components();
+
+ Assert (first_comp < fe.n_components(), ExcInternalError());
+ Assert (end_comp <= fe.n_components(), ExcInternalError());
+
+ // now check whether any of
+ // the components in between
+ // is set
+ for (unsigned int c=first_comp; c<end_comp; ++c)
+ if (component_select[c] == true)
+ {
+ local_selected_dofs[i] = true;
+ break;
+ }
+ }
// then loop over all cells and do
// the work
template<int dim>
void
-DoFTools::extract_level_dofs(const unsigned int level,
- const MGDoFHandler<dim> &dof,
- const std::vector<bool> &component_select,
- std::vector<bool> &selected_dofs)
+DoFTools::extract_level_dofs(
+ const unsigned int level,
+ const MGDoFHandler<dim> &dof,
+ const std::vector<bool> &component_select,
+ std::vector<bool> &selected_dofs,
+ bool count_by_blocks)
{
const FiniteElement<dim>& fe = dof.get_fe();
- Assert(component_select.size() == fe.n_components(),
- ExcDimensionMismatch(component_select.size(), fe.n_components()));
+
+ if (count_by_blocks)
+ {
+ Assert(component_select.size() == fe.n_blocks(),
+ ExcDimensionMismatch(component_select.size(), fe.n_blocks()));
+ }
+ else
+ {
+ Assert(component_select.size() == fe.n_components(),
+ ExcDimensionMismatch(component_select.size(), fe.n_components()));
+ }
+
Assert(selected_dofs.size() == dof.n_dofs(level),
ExcDimensionMismatch(selected_dofs.size(), dof.n_dofs(level)));
// something interesting or not
std::vector<bool> local_selected_dofs (fe.dofs_per_cell, false);
for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
- if (fe.is_primitive(i))
+ if (count_by_blocks)
local_selected_dofs[i]
- = component_select[fe.system_to_component_index(i).first];
+ = component_select[fe.system_to_block_index(i).first];
else
- // if this shape function is
- // not primitive, then we have
- // to work harder. we have to
- // find out whether _any_ of
- // the vector components of
- // this element is selected or
- // not
- //
- // to do so, get the first and
- // last vector components of
- // the base element to which
- // the local dof with index i
- // belongs
- {
- unsigned int first_comp = 0;
- const unsigned int this_base = fe.system_to_base_index(i).first.first;
- const unsigned int this_multiplicity
- = fe.system_to_base_index(i).first.second;
-
- for (unsigned int b=0; b<this_base; ++b)
- first_comp += fe.base_element(b).n_components() *
- fe.element_multiplicity(b);
- for (unsigned int m=0; m<this_multiplicity; ++m)
- first_comp += fe.base_element(this_base).n_components();
- const unsigned int end_comp = first_comp +
- fe.base_element(this_base).n_components();
-
- Assert (first_comp < fe.n_components(), ExcInternalError());
- Assert (end_comp <= fe.n_components(), ExcInternalError());
-
- // now check whether any of
- // the components in between
- // is set
- for (unsigned int c=first_comp; c<end_comp; ++c)
- if (component_select[c] == true)
- {
- local_selected_dofs[i] = true;
- break;
- };
- };
-
+ if (fe.is_primitive(i))
+ local_selected_dofs[i]
+ = component_select[fe.system_to_component_index(i).first];
+ else
+ // if this shape function is
+ // not primitive, then we have
+ // to work harder. we have to
+ // find out whether _any_ of
+ // the vector components of
+ // this element is selected or
+ // not
+ //
+ // to do so, get the first and
+ // last vector components of
+ // the base element to which
+ // the local dof with index i
+ // belongs
+ {
+ unsigned int first_comp = 0;
+ const unsigned int this_base = fe.system_to_base_index(i).first.first;
+ const unsigned int this_multiplicity
+ = fe.system_to_base_index(i).first.second;
+
+ for (unsigned int b=0; b<this_base; ++b)
+ first_comp += fe.base_element(b).n_components() *
+ fe.element_multiplicity(b);
+ for (unsigned int m=0; m<this_multiplicity; ++m)
+ first_comp += fe.base_element(this_base).n_components();
+ const unsigned int end_comp = first_comp +
+ fe.base_element(this_base).n_components();
+
+ Assert (first_comp < fe.n_components(), ExcInternalError());
+ Assert (end_comp <= fe.n_components(), ExcInternalError());
+
+ // now check whether any of
+ // the components in between
+ // is set
+ for (unsigned int c=first_comp; c<end_comp; ++c)
+ if (component_select[c] == true)
+ {
+ local_selected_dofs[i] = true;
+ break;
+ }
+ }
+
// then loop over all cells and do
// work
std::vector<unsigned int> indices(fe.dofs_per_cell);
{
void (*fun_ptr) (const DoFHandler<dim> &,
const std::vector<bool> &,
- std::vector<bool> &)
+ std::vector<bool> &,
+ bool)
= &DoFTools::template extract_dofs<dim>;
component_select[i][i] = true;
threads += Threads::spawn (fun_ptr)(dof_handler, component_select[i],
- dofs_in_component[i]);
+ dofs_in_component[i], false);
};
threads.join_all ();
template void DoFTools::extract_dofs<deal_II_dimension>
-(const DoFHandler<deal_II_dimension>& dof,
- const std::vector<bool>& component_select,
- std::vector<bool>& selected_dofs);
+(const DoFHandler<deal_II_dimension>&,
+ const std::vector<bool>&, std::vector<bool>&, bool);
template void DoFTools::extract_level_dofs<deal_II_dimension>
-(const unsigned int level,
- const MGDoFHandler<deal_II_dimension>& dof,
- const std::vector<bool>& component_select,
- std::vector<bool>& selected_dofs);
+(const unsigned int level, const MGDoFHandler<deal_II_dimension>&,
+ const std::vector<bool>&, std::vector<bool>&, bool);
#if deal_II_dimension != 1
template