]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
VectorTools::interpolation_to_different_mesh implemented
authormaier <maier@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 8 Feb 2012 14:40:45 +0000 (14:40 +0000)
committermaier <maier@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 8 Feb 2012 14:40:45 +0000 (14:40 +0000)
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

deal.II/doc/news/changes.h
deal.II/include/deal.II/numerics/vectors.h
deal.II/include/deal.II/numerics/vectors.templates.h
deal.II/source/numerics/vectors.inst.in

index 06441b12dc7ea617a6d39667817865b72df847b7..059568060c416b22eb9273acbbda8a83ead3878d 100644 (file)
@@ -87,6 +87,12 @@ enabled due to a missing include file in file
 <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)
index 0bf46f634bb8e6b3b5c6579fbfea12d0783b10ee..d34c18b8884e9eaf22044856fd723f07f8a1cf44 100644 (file)
@@ -36,6 +36,7 @@ template <typename number> class Vector;
 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;
@@ -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 <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.
index 0bb918728339b010e7913ce977b3a287e47168f3..8cdbfbfb92b01897d1dd036180b436ce7c9f3b00 100644 (file)
@@ -34,6 +34,7 @@
 #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>
@@ -421,6 +422,112 @@ namespace VectorTools
 
 
 
+  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,
index 0f98a636189d9c09330084eeb29e7738f5ae8366..ad50778264dc85bf867af74038b1fe5c168e7293 100644 (file)
@@ -299,6 +299,7 @@ for (VEC : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS)
      const unsigned int);
 #endif
 
+
   template
     void project
     (const Mapping<deal_II_dimension>      &,
@@ -838,3 +839,56 @@ void create_right_hand_side<deal_II_dimension,deal_II_dimension+1>
   \}
 }
 
+
+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>                                     &);
+
+  \}
+}
+
+

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.