(Timo Heister, 2014/01/17)
</li>
- <li> Fixed: DoFCellAccessor::set_dof_values_by_interpolation and
+ <li> Fixed/new: DoFCellAccessor::set_dof_values_by_interpolation and
DoFCellAccessor::get_interpolated_dof_values could previously be
called for hp::DoFHandler objects on cells that are non-active. This
makes no sense since these cells have no associated finite element
space. Doing so now raises an exception.
<br>
- (Wolfgang Bangerth, 2014/01/17)
+ However, there are legitimate cases where one may want to interpolate
+ from children to a parent's finite element space. Since in the hp
+ case no finite element space is naturally associated with an inactive
+ cell, it is now possible to pass an explicit finite element index
+ argument to the function specifying which element of an hp::FECollection
+ object describes the space onto which you want to interpolate.
+ <br>
+ (Wolfgang Bangerth, 2014/01/18)
</li>
<li> Fixed: The methods IndexSet::do_compress() and
/**
* Return the interpolation of the given finite element function to
* the present cell. In the simplest case, the cell is a terminal
- * one, i.e. has no children; then, the returned value is the vector
- * of nodal values on that cell. You could then as well get the
+ * one, i.e., it has no children; then, the returned value is the vector
+ * of nodal values on that cell. You could as well get the
* desired values through the @p get_dof_values function. In the
* other case, when the cell has children, we use the restriction
* matrices provided by the finite element class to compute the
* interpolation from the children to the present cell.
*
- * It is assumed that both vectors already have the right size
+ * If the cell is part of a hp::DoFHandler object, cells only have an
+ * associated finite element space if they are active. However, this
+ * function is supposed to also provide information on inactive cells
+ * with children. Consequently, it carries a third argument that can be
+ * used in the hp context that denotes the finite element space we are
+ * supposed to interpolate onto. If the cell is active, this function
+ * then obtains the finite element function from the <code>values</code>
+ * vector on this cell and interpolates it onto the space described
+ * by the <code>fe_index</code>th element of the hp::FECollection associated
+ * with the hp::DoFHandler of which this cell is a part of. If the cell
+ * is not active, then we first perform this interpolation on all of its
+ * terminal children and then interpolate this function down to the
+ * cell requested keeping the function space the same.
+ *
+ * It is assumed that both input vectors already have the right size
* beforehand.
*
- * Unlike the get_dof_values() function, this function works on
- * cells rather than to lines, quads, and hexes, since interpolation
+ * @note Unlike the get_dof_values() function, this function is only available
+ * on cells, rather than on lines, quads, and hexes, since interpolation
* is presently only provided for cells by the finite element
* classes.
*/
template <class InputVector, typename number>
void get_interpolated_dof_values (const InputVector &values,
- Vector<number> &interpolated_values) const;
+ Vector<number> &interpolated_values,
+ const unsigned int fe_index = DH::default_fe_index) const;
/**
- * This, again, is the counterpart to get_interpolated_dof_values():
+ * This function is the counterpart to get_interpolated_dof_values():
* you specify the dof values on a cell and these are interpolated
* to the children of the present cell and set on the terminal
* cells.
*/
template <class OutputVector, typename number>
void set_dof_values_by_interpolation (const Vector<number> &local_values,
- OutputVector &values) const;
+ OutputVector &values,
+ const unsigned int fe_index = DH::default_fe_index) const;
/**
* Distribute a local (cell based) vector to a global one by mapping
void
DoFCellAccessor<DH,lda>::
get_interpolated_dof_values (const InputVector &values,
- Vector<number> &interpolated_values) const
+ Vector<number> &interpolated_values,
+ const unsigned int fe_index) const
{
if (!this->has_children())
// if this cell has no children: simply return the exact values on this
- // cell
- this->get_dof_values (values, interpolated_values);
+ // cell unless the finite element we need to interpolate to is different than
+ // the one we have on the current cell
+ {
+ if ((dynamic_cast<DoFHandler<DH::dimension,DH::space_dimension>*>
+ (this->dof_handler)
+ != 0)
+ ||
+ (fe_index == this->active_fe_index())
+ ||
+ (fe_index == DH::default_fe_index))
+ this->get_dof_values (values, interpolated_values);
+ else
+ {
+ // well, here we need to first get the values from the current
+ // cell and then interpolate it to the element requested. this
+ // can clearly only happen for hp::DoFHandler objects
+ Vector<number> tmp (this->get_fe().dofs_per_cell);
+ this->get_dof_values (values, tmp);
+
+ FullMatrix<double> interpolation (this->dof_handler->get_fe()[fe_index].dofs_per_cell,
+ this->get_fe().dofs_per_cell);
+ this->dof_handler->get_fe()[fe_index].get_interpolation_matrix (this->get_fe(),
+ interpolation);
+ interpolation.vmult (interpolated_values, tmp);
+ }
+ }
else
// otherwise obtain them from the children
{
// context, we can simply assume that the FE space to which we
// want to interpolate is the same as for all elements in the
// mesh). consequently, we cannot interpolate from children's FE
- // space to this cell's (unknown) FE space
+ // space to this cell's (unknown) FE space unless an explicit
+ // fe_index is given
Assert ((dynamic_cast<DoFHandler<DH::dimension,DH::space_dimension>*>
(this->dof_handler)
- != 0),
+ != 0)
+ ||
+ (fe_index != DH::default_fe_index),
ExcMessage ("You cannot call this function on non-active cells "
- "of hp::DoFHandler objects because they do not have "
+ "of hp::DoFHandler objects unless you provide an explicit "
+ "finite element index because they do not have naturally "
"associated finite element spaces associated: degrees "
"of freedom are only distributed on active cells for which "
"the active_fe_index has been set."));
- const FiniteElement<dim,spacedim> &fe = this->get_fe();
+ const FiniteElement<dim,spacedim> &fe = this->get_dof_handler().get_fe()[fe_index];
const unsigned int dofs_per_cell = fe.dofs_per_cell;
Assert (this->dof_handler != 0,
for (unsigned int child=0; child<this->n_children(); ++child)
{
// get the values from the present child, if necessary by
- // interpolation itself
+ // interpolation itself either from its own children or
+ // by interpolating from the finite element on an active
+ // child to the finite element space requested here
this->child(child)->get_interpolated_dof_values (values,
- tmp1);
+ tmp1,
+ fe_index);
// interpolate these to the mother cell
fe.get_restriction_matrix(child, this->refinement_case()).vmult (tmp2, tmp1);
template
void
DoFCellAccessor<DoFHandler<deal_II_dimension>, lda>::get_interpolated_dof_values
- (const VEC&, Vector<SCALAR>&) const;
+ (const VEC&, Vector<SCALAR>&,
+ const unsigned int fe_index) const;
#if deal_II_dimension != 3
template
void
DoFCellAccessor<DoFHandler<deal_II_dimension,deal_II_dimension+1>, lda>::get_interpolated_dof_values
- (const VEC&, Vector<SCALAR>&) const;
+ (const VEC&, Vector<SCALAR>&,
+ const unsigned int fe_index) const;
#endif
template
void
DoFCellAccessor<DoFHandler<1,3>, lda>::get_interpolated_dof_values
- (const VEC&, Vector<SCALAR>&) const;
+ (const VEC&, Vector<SCALAR>&,
+ const unsigned int fe_index) const;
#endif
template
void
DoFCellAccessor<hp::DoFHandler<deal_II_dimension>, lda>::get_interpolated_dof_values
- (const VEC&, Vector<SCALAR>&) const;
+ (const VEC&, Vector<SCALAR>&,
+ const unsigned int fe_index) const;
#if deal_II_dimension != 3
template
void
DoFCellAccessor<hp::DoFHandler<deal_II_dimension,deal_II_dimension+1>, lda>::get_interpolated_dof_values
- (const VEC&, Vector<SCALAR>&) const;
+ (const VEC&, Vector<SCALAR>&,
+ const unsigned int fe_index) const;
#endif
template
void
DoFCellAccessor<hp::DoFHandler<1,3>, lda>::get_interpolated_dof_values
- (const VEC&, Vector<SCALAR>&) const;
+ (const VEC&, Vector<SCALAR>&,
+ const unsigned int fe_index) const;
#endif
}
void
DoFCellAccessor<DH,lda>::
set_dof_values_by_interpolation (const Vector<number> &local_values,
- OutputVector &values) const
+ OutputVector &values,
+ const unsigned int fe_index) const
{
if (!this->has_children())
// if this cell has no children: simply set the values on this cell
// context, we can simply assume that the FE space to from which we
// want to interpolate is the same as for all elements in the
// mesh). consequently, we cannot interpolate to children's FE
- // space from this cell's (unknown) FE space
+ // space from this cell's (unknown) FE space unless an explicit
+ // fe_index is given
Assert ((dynamic_cast<DoFHandler<DH::dimension,DH::space_dimension>*>
- (this->dof_handler)
- != 0),
- ExcMessage ("You cannot call this function on non-active cells "
- "of hp::DoFHandler objects because they do not have "
- "associated finite element spaces associated: degrees "
- "of freedom are only distributed on active cells for which "
- "the active_fe_index has been set."));
+ (this->dof_handler)
+ != 0)
+ ||
+ (fe_index != DH::default_fe_index),
+ ExcMessage ("You cannot call this function on non-active cells "
+ "of hp::DoFHandler objects unless you provide an explicit "
+ "finite element index because they do not have naturally "
+ "associated finite element spaces associated: degrees "
+ "of freedom are only distributed on active cells for which "
+ "the active_fe_index has been set."));
const unsigned int dofs_per_cell = this->get_fe().dofs_per_cell;
template
void
DoFCellAccessor<DoFHandler<deal_II_dimension>, lda>::set_dof_values_by_interpolation
- (const Vector<SCALAR>&, VEC&) const;
+ (const Vector<SCALAR>&, VEC&,
+ const unsigned int fe_index) const;
#if deal_II_dimension != 3
template
void
DoFCellAccessor<DoFHandler<deal_II_dimension,deal_II_dimension+1>, lda>::set_dof_values_by_interpolation
- (const Vector<SCALAR>&, VEC&) const;
+ (const Vector<SCALAR>&, VEC&,
+ const unsigned int fe_index) const;
#endif
template
void
DoFCellAccessor<DoFHandler<1,3>, lda>::set_dof_values_by_interpolation
- (const Vector<SCALAR>&, VEC&) const;
+ (const Vector<SCALAR>&, VEC&,
+ const unsigned int fe_index) const;
#endif
template
void
DoFCellAccessor<hp::DoFHandler<deal_II_dimension>, lda>::set_dof_values_by_interpolation
- (const Vector<SCALAR>&, VEC&) const;
+ (const Vector<SCALAR>&, VEC&,
+ const unsigned int fe_index) const;
#if deal_II_dimension != 3
template
void
DoFCellAccessor<hp::DoFHandler<deal_II_dimension,deal_II_dimension+1>, lda>::set_dof_values_by_interpolation
- (const Vector<SCALAR>&, VEC&) const;
+ (const Vector<SCALAR>&, VEC&,
+ const unsigned int fe_index) const;
#endif
template
void
DoFCellAccessor<hp::DoFHandler<1,3>, lda>::set_dof_values_by_interpolation
- (const Vector<SCALAR>&, VEC&) const;
+ (const Vector<SCALAR>&, VEC&,
+ const unsigned int fe_index) const;
#endif
}