<h3>Specific improvements</h3>
<ol>
+
+<li> New: VectorTools::interpolation_to_different_mesh implemented which interpolates between
+ DoFHandlers with different triangulations based on a common coarse grid.
+<br>
+(Matthias Maier, 2012/02/08)
+
<li> Improved: DoFTools::map_dofs_to_support_points() now also works within the hp framework.
<br>
(Christian Goll, 2012/02/02)
template <typename number> class FullMatrix;
template <int dim, int spacedim> class Mapping;
template <int dim, int spacedim> class DoFHandler;
+template <typename gridtype> class InterGridMap;
namespace hp
{
template <int dim, int spacedim> class DoFHandler;
const InVector &data_1,
OutVector &data_2);
+ /**
+ * Gives the interpolation of a
+ * @p dof1-function @p u1 to a
+ * @p dof2-function @p u2, where @p
+ * dof1 and @p dof2 represent
+ * different triangulations with a
+ * common coarse grid.
+ *
+ * dof1 and dof2 need to have the
+ * same finite element
+ * discretization.
+ *
+ * Note that for continuous
+ * elements on grids with hanging
+ * nodes (i.e. locally refined
+ * grids) this function does not
+ * give the expected output.
+ * Indeed, the resulting output
+ * vector does not necessarily
+ * respect continuity
+ * requirements at hanging nodes,
+ * due to local cellwise
+ * interpolation.
+ *
+ * For this case (continuous
+ * elements on grids with hanging
+ * nodes), please use the
+ * interpolate_to_different_mesh
+ * function with an additional
+ * ConstraintMatrix argument,
+ * see below, or make the field
+ * conforming yourself by calling the
+ * @p ConstraintsMatrix::distribute
+ * function of your hanging node
+ * constraints object.
+ */
+ template <int dim, int spacedim,
+ template <int,int> class DH,
+ class VECTOR>
+ void
+ interpolate_to_different_mesh (const DH<dim, spacedim> &dof1,
+ const VECTOR &u1,
+ const DH<dim, spacedim> &dof2,
+ VECTOR &u2);
+
+ /**
+ * Gives the interpolation of a
+ * @p dof1-function @p u1 to a
+ * @p dof2-function @p u2, where @p
+ * dof1 and @p dof2 represent
+ * different triangulations with a
+ * common coarse grid.
+ *
+ * dof1 and dof2 need to have the
+ * same finite element
+ * discretization.
+ *
+ * @p constraints is a hanging node
+ * constraints object corresponding
+ * to @p dof2. This object is
+ * particularly important when
+ * interpolating onto continuous
+ * elements on grids with hanging
+ * nodes (locally refined grids):
+ * Without it - due to cellwise
+ * interpolation - the resulting
+ * output vector does not necessarily
+ * respect continuity requirements
+ * at hanging nodes.
+ */
+ template <int dim, int spacedim,
+ template <int,int> class DH,
+ class VECTOR>
+ void
+ interpolate_to_different_mesh (const DH<dim, spacedim> &dof1,
+ const VECTOR &u1,
+ const DH<dim, spacedim> &dof2,
+ const ConstraintMatrix &constraints,
+ VECTOR &u2);
+
+
+ /**
+ * The same function as above, but
+ * takes an InterGridMap object
+ * directly as a parameter. Useful
+ * for interpolating several vectors
+ * at the same time.
+ *
+ * @p intergridmap
+ * has to be initialized via
+ * InterGridMap::make_mapping pointing
+ * from a source DoFHandler to a
+ * destination DoFHandler.
+ */
+ template <int dim, int spacedim,
+ template <int,int> class DH,
+ class VECTOR>
+ void
+ interpolate_to_different_mesh (const InterGridMap<DH<dim, spacedim> > &intergridmap,
+ const VECTOR &u1,
+ const ConstraintMatrix &constraints,
+ VECTOR &u2);
+
/**
* Compute the projection of
* @p function to the finite element space.
#include <deal.II/grid/tria_accessor.h>
#include <deal.II/grid/tria_boundary.h>
#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/intergrid_map.h>
#include <deal.II/dofs/dof_accessor.h>
#include <deal.II/dofs/dof_handler.h>
#include <deal.II/dofs/dof_tools.h>
+ template <int dim, int spacedim,
+ template <int,int> class DH,
+ class Vector>
+ void
+ interpolate_to_different_mesh (const DH<dim, spacedim> &dof1,
+ const Vector &u1,
+ const DH<dim, spacedim> &dof2,
+ Vector &u2)
+ {
+ Assert(GridTools::have_same_coarse_mesh(dof1, dof2),
+ ExcMessage ("The two containers must represent triangulations that "
+ "have the same coarse meshes"));
+
+ InterGridMap<DH<dim, spacedim> > intergridmap;
+ intergridmap.make_mapping(dof1, dof2);
+
+ ConstraintMatrix dummy;
+ dummy.close();
+
+ interpolate_to_different_mesh(intergridmap, u1, dummy, u2);
+ }
+
+
+
+ template <int dim, int spacedim,
+ template <int,int> class DH,
+ class Vector>
+ void
+ interpolate_to_different_mesh (const DH<dim, spacedim> &dof1,
+ const Vector &u1,
+ const DH<dim, spacedim> &dof2,
+ const ConstraintMatrix &constraints,
+ Vector &u2)
+ {
+ Assert(GridTools::have_same_coarse_mesh(dof1, dof2),
+ ExcMessage ("The two containers must represent triangulations that "
+ "have the same coarse meshes"));
+
+ InterGridMap<DH<dim, spacedim> > intergridmap;
+ intergridmap.make_mapping(dof1, dof2);
+
+ interpolate_to_different_mesh(intergridmap, u1, constraints, u2);
+ }
+
+
+
+ template <int dim, int spacedim,
+ template <int,int> class DH,
+ class Vector>
+ void
+ interpolate_to_different_mesh (const InterGridMap<DH<dim, spacedim> > &intergridmap,
+ const Vector &u1,
+ const ConstraintMatrix &constraints,
+ Vector &u2)
+ {
+ const DH<dim, spacedim> &dof1 = intergridmap.get_source_grid();
+ const DH<dim, spacedim> &dof2 = intergridmap.get_destination_grid();
+
+ Assert(u1.size()==dof1.n_dofs(),
+ ExcDimensionMismatch(u1.size(), dof1.n_dofs()));
+ Assert(u2.size()==dof2.n_dofs(),
+ ExcDimensionMismatch(u2.size(),dof2.n_dofs()));
+
+ Vector cache;
+
+ // Looping over the finest common
+ // mesh, this means that source and
+ // destination cells have to be on the
+ // same level and at least one has to
+ // be active.
+ //
+ // Therefor, loop over all cells
+ // (active and inactive) of the source
+ // grid ..
+ typename DH<dim,spacedim>::cell_iterator cell1 = dof1.begin();
+ const typename DH<dim,spacedim>::cell_iterator endc1 = dof1.end();
+
+ for (; cell1 != endc1; ++cell1) {
+ const typename DH<dim,spacedim>::cell_iterator cell2 = intergridmap[cell1];
+
+ // .. and skip if source and destination
+ // cells are not on the same level ..
+ if (cell1->level() != cell2->level())
+ continue;
+ // .. or none of them is active.
+ if (!cell1->active() && !cell2->active())
+ continue;
+
+ Assert(cell1->get_fe().get_name() ==
+ cell2->get_fe().get_name(),
+ ExcMessage ("Source and destination cells need to use the same finite element"));
+
+ cache.reinit(dof1.get_fe().dofs_per_cell);
+
+ // Get and set the corresponding
+ // dof_values by interpolation.
+ cell1->get_interpolated_dof_values(u1, cache);
+ cell2->set_dof_values_by_interpolation(cache, u2);
+ }
+
+ // Apply hanging node constraints.
+ constraints.distribute(u2);
+ }
+
+
+
template <int dim, class VECTOR, int spacedim>
void project (const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
const unsigned int);
#endif
+
template
void project
(const Mapping<deal_II_dimension> &,
\}
}
+
+for (deal_II_dimension : DIMENSIONS)
+{
+ namespace VectorTools \{
+
+ template
+ void interpolate_to_different_mesh
+ (const DoFHandler<deal_II_dimension> &,
+ const Vector<float> &,
+ const DoFHandler<deal_II_dimension> &,
+ Vector<float> &);
+
+ template
+ void interpolate_to_different_mesh
+ (const DoFHandler<deal_II_dimension> &,
+ const Vector<float> &,
+ const DoFHandler<deal_II_dimension> &,
+ const ConstraintMatrix &,
+ Vector<float> &);
+
+ template
+ void interpolate_to_different_mesh
+ (const InterGridMap<DoFHandler<deal_II_dimension> > &,
+ const Vector<float> &,
+ const ConstraintMatrix &,
+ Vector<float> &);
+
+ template
+ void interpolate_to_different_mesh
+ (const DoFHandler<deal_II_dimension> &,
+ const Vector<double> &,
+ const DoFHandler<deal_II_dimension> &,
+ Vector<double> &);
+
+ template
+ void interpolate_to_different_mesh
+ (const DoFHandler<deal_II_dimension> &,
+ const Vector<double> &,
+ const DoFHandler<deal_II_dimension> &,
+ const ConstraintMatrix &,
+ Vector<double> &);
+
+ template
+ void interpolate_to_different_mesh
+ (const InterGridMap<DoFHandler<deal_II_dimension> > &,
+ const Vector<double> &,
+ const ConstraintMatrix &,
+ Vector<double> &);
+
+ \}
+}
+
+