* Note that the resulting field
* does not satisfy continuity
* requirements of the given
- * finite elements. You have to
- * apply
- * @ref{ConstraintMatrix}@p{::distribute}
- * with the hanging node
- * constraints of the second DoF
- * handler object to make the
- * output continuous again.
+ * finite elements.
+ *
+ * When you use continuous
+ * elements on grids with hanging
+ * nodes, please use the
+ * @p{extrapolate} function with
+ * an additional
+ * @p{ConstraintMatrix} argument,
+ * see below.
+ */
+ template <int dim, typename number>
+ static void extrapolate(const DoFHandler<dim> &dof1,
+ const Vector<number> &z1,
+ const DoFHandler<dim> &dof2,
+ Vector<number> &z2);
+
+ /**
+ * Gives the patchwise
+ * extrapolation of a @p{dof1}
+ * function @p{z1} to a @p{dof2}
+ * function @p{z2}. @p{dof1} and
+ * @p{dof2} need to be @ref{DoFHandler}
+ * based on the same triangulation.
+ * @p{constraints} is a hanging
+ * node constraints object
+ * corresponding to
+ * @p{dof2}. This object is
+ * particular important when
+ * interpolating onto continuous
+ * elements on grids with hanging
+ * nodes (locally refined grids).
+ *
+ * This function is interesting
+ * for e.g. extrapolating
+ * patchwise a piecewise linear
+ * solution to a piecewise
+ * quadratic solution.
*/
template <int dim, typename number>
static void extrapolate(const DoFHandler<dim> &dof1,
const Vector<number> &z1,
const DoFHandler<dim> &dof2,
+ const ConstraintMatrix &constraints,
Vector<number> &z2);
/**
const Vector<number> &u1,
const DoFHandler<dim> &dof2,
Vector<number> &u2)
+{
+ ConstraintMatrix dummy;
+ dummy.close();
+ extrapolate(dof1, u1, dof2, dummy, u2);
+}
+
+
+
+template <int dim, typename number>
+void FETools::extrapolate(const DoFHandler<dim> &dof1,
+ const Vector<number> &u1,
+ const DoFHandler<dim> &dof2,
+ const ConstraintMatrix &constraints,
+ Vector<number> &u2)
{
Assert(dof1.get_fe().n_components() == dof2.get_fe().n_components(),
ExcDimensionMismatch(dof1.get_fe().n_components(), dof2.get_fe().n_components()));
Assert(u2.size()==dof2.n_dofs(), ExcDimensionMismatch(u2.size(), dof2.n_dofs()));
Vector<number> u3(dof2.n_dofs());
- interpolate(dof1, u1, dof2, u3);
+ interpolate(dof1, u1, dof2, constraints, u3);
const unsigned int dofs_per_cell = dof2.get_fe().dofs_per_cell;
const Triangulation<dim> &tria=dof1.get_tria();
}
}
}
+
+ // Apply hanging node constraints.
+ constraints.distribute(u2);
}
<h3>deal.II</h3>
<ol>
+ <li> <p> New: The existing <code
+ class="member">FETools::extrapolate</code> functions does not
+ take into account hanging node constraints. Therefore, it works
+ for continuous elements on globally refined grids and on
+ discontinuous elements, only. Now, there is a new <code
+ class="member">FETools::extrapolate</code> function with an
+ additional <code class="class">ConstraintMatrix</code> argument
+ for the hanging node constraints. This new function works for
+ continuous as well as for discontinuous elements. But, the old
+ function is still supported.
+ <br>
+ (RH 2002/06/17)
+ </p>
+
<li> <p> New: There are now <code
class="member">Triangulation::save/load_user_pointers(vector<void
*> &)</code> functions similar to the respective <code