]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
New extrapolate function taking into account hanging node constraints.
authorhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 17 Jun 2002 14:21:54 +0000 (14:21 +0000)
committerhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 17 Jun 2002 14:21:54 +0000 (14:21 +0000)
git-svn-id: https://svn.dealii.org/trunk@6151 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_tools.h
deal.II/deal.II/source/fe/fe_tools.cc
deal.II/doc/news/2002/c-3-4.html

index c4cfebc61f8b448ce17ee13860fcf1facb25d8cc..d3b0b534c119c17c82381681c6c54f52e418478a 100644 (file)
@@ -341,18 +341,49 @@ class FETools
                                      * 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);    
 
                                     /**
index 5ad56f1cdb8f7b2ea87b8bf71c276e2aca9a1ea6..3bbfb5b61a1672669882bff75f9074f9a76f2b39 100644 (file)
@@ -375,6 +375,20 @@ void FETools::extrapolate(const DoFHandler<dim> &dof1,
                          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()));
@@ -383,7 +397,7 @@ void FETools::extrapolate(const DoFHandler<dim> &dof1,
   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();
@@ -423,6 +437,9 @@ void FETools::extrapolate(const DoFHandler<dim> &dof1,
              }
          }
     }
+
+                                  // Apply hanging node constraints.
+  constraints.distribute(u2);
 }
 
 
index 04447c0734e5d8fe2831ca92f65f697a448dfd8b..c1d0f0c879132145d3df2d44d5cdccb6e684ed33 100644 (file)
@@ -105,6 +105,20 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
 <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&lt;void
        *&gt; &)</code> functions similar to the respective <code

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.