From: maier Date: Wed, 8 Feb 2012 14:40:45 +0000 (+0000) Subject: VectorTools::interpolation_to_different_mesh implemented X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=c5a220254c2056e217887d9ce8820fe52e97c0cf;p=dealii-svn.git VectorTools::interpolation_to_different_mesh implemented VectorTools::interpolation_to_different_mesh gives the interpolation of a dof1-function u1 to a dof2-function u2, where dof1 and dof2 represent different triangulations with a common coarse grid. dof1 and dof2 need to have the same finite element discretization. The mapping between the two spaces is done via the help of the InterGridMap<> class. git-svn-id: https://svn.dealii.org/trunk@25013 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 06441b12dc..059568060c 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -87,6 +87,12 @@ enabled due to a missing include file in file

Specific improvements

    + +
  1. New: VectorTools::interpolation_to_different_mesh implemented which interpolates between + DoFHandlers with different triangulations based on a common coarse grid. +
    +(Matthias Maier, 2012/02/08) +
  2. Improved: DoFTools::map_dofs_to_support_points() now also works within the hp framework.
    (Christian Goll, 2012/02/02) diff --git a/deal.II/include/deal.II/numerics/vectors.h b/deal.II/include/deal.II/numerics/vectors.h index 0bf46f634b..d34c18b888 100644 --- a/deal.II/include/deal.II/numerics/vectors.h +++ b/deal.II/include/deal.II/numerics/vectors.h @@ -36,6 +36,7 @@ template class Vector; template class FullMatrix; template class Mapping; template class DoFHandler; +template class InterGridMap; namespace hp { template class DoFHandler; @@ -495,6 +496,109 @@ namespace VectorTools 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 class DH, + class VECTOR> + void + interpolate_to_different_mesh (const DH &dof1, + const VECTOR &u1, + const DH &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 class DH, + class VECTOR> + void + interpolate_to_different_mesh (const DH &dof1, + const VECTOR &u1, + const DH &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 class DH, + class VECTOR> + void + interpolate_to_different_mesh (const InterGridMap > &intergridmap, + const VECTOR &u1, + const ConstraintMatrix &constraints, + VECTOR &u2); + /** * Compute the projection of * @p function to the finite element space. diff --git a/deal.II/include/deal.II/numerics/vectors.templates.h b/deal.II/include/deal.II/numerics/vectors.templates.h index 0bb9187283..8cdbfbfb92 100644 --- a/deal.II/include/deal.II/numerics/vectors.templates.h +++ b/deal.II/include/deal.II/numerics/vectors.templates.h @@ -34,6 +34,7 @@ #include #include #include +#include #include #include #include @@ -421,6 +422,112 @@ namespace VectorTools + template class DH, + class Vector> + void + interpolate_to_different_mesh (const DH &dof1, + const Vector &u1, + const DH &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 > intergridmap; + intergridmap.make_mapping(dof1, dof2); + + ConstraintMatrix dummy; + dummy.close(); + + interpolate_to_different_mesh(intergridmap, u1, dummy, u2); + } + + + + template class DH, + class Vector> + void + interpolate_to_different_mesh (const DH &dof1, + const Vector &u1, + const DH &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 > intergridmap; + intergridmap.make_mapping(dof1, dof2); + + interpolate_to_different_mesh(intergridmap, u1, constraints, u2); + } + + + + template class DH, + class Vector> + void + interpolate_to_different_mesh (const InterGridMap > &intergridmap, + const Vector &u1, + const ConstraintMatrix &constraints, + Vector &u2) + { + const DH &dof1 = intergridmap.get_source_grid(); + const DH &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::cell_iterator cell1 = dof1.begin(); + const typename DH::cell_iterator endc1 = dof1.end(); + + for (; cell1 != endc1; ++cell1) { + const typename DH::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 void project (const Mapping &mapping, const DoFHandler &dof, diff --git a/deal.II/source/numerics/vectors.inst.in b/deal.II/source/numerics/vectors.inst.in index 0f98a63618..ad50778264 100644 --- a/deal.II/source/numerics/vectors.inst.in +++ b/deal.II/source/numerics/vectors.inst.in @@ -299,6 +299,7 @@ for (VEC : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS) const unsigned int); #endif + template void project (const Mapping &, @@ -838,3 +839,56 @@ void create_right_hand_side \} } + +for (deal_II_dimension : DIMENSIONS) +{ + namespace VectorTools \{ + + template + void interpolate_to_different_mesh + (const DoFHandler &, + const Vector &, + const DoFHandler &, + Vector &); + + template + void interpolate_to_different_mesh + (const DoFHandler &, + const Vector &, + const DoFHandler &, + const ConstraintMatrix &, + Vector &); + + template + void interpolate_to_different_mesh + (const InterGridMap > &, + const Vector &, + const ConstraintMatrix &, + Vector &); + + template + void interpolate_to_different_mesh + (const DoFHandler &, + const Vector &, + const DoFHandler &, + Vector &); + + template + void interpolate_to_different_mesh + (const DoFHandler &, + const Vector &, + const DoFHandler &, + const ConstraintMatrix &, + Vector &); + + template + void interpolate_to_different_mesh + (const InterGridMap > &, + const Vector &, + const ConstraintMatrix &, + Vector &); + + \} +} + +