]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement VectorTools::project also for 1d.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 9 Aug 2006 03:28:13 +0000 (03:28 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 9 Aug 2006 03:28:13 +0000 (03:28 +0000)
git-svn-id: https://svn.dealii.org/trunk@13620 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/vectors.h
deal.II/deal.II/include/numerics/vectors.templates.h
deal.II/doc/news/changes.html

index b50bf9c89bb2ca72633b1334164533d55955efee..ca7561820ab3eb2eaa9b197f437ac53e41a671fb 100644 (file)
@@ -496,6 +496,13 @@ class VectorTools
                                      *
                                      * See the general documentation of this
                                      * class for further information.
+                                     * 
+                                     * In 1d, the default value of
+                                     * the boundary quadrature
+                                     * formula is an invalid object
+                                     * since integration on the
+                                     * boundary doesn't happen in
+                                     * 1d.
                                      */
     template <int dim, class VECTOR>
     static void project (const Mapping<dim>       &mapping,
@@ -505,28 +512,9 @@ class VectorTools
                         const Function<dim>      &function,
                         VECTOR                   &vec,
                         const bool                enforce_zero_boundary = false,
-                        const Quadrature<dim-1>  &q_boundary = QGauss<dim-1>(2),
-                        const bool                project_to_boundary_first = false);
-
-                                    /**
-                                     * Declaration of specialization
-                                     * of the previous function for
-                                     * 1d. At present, it is not
-                                     * implemented.
-                                     *
-                                     * The default value of the boundary
-                                     * quadrature formula is an invalid
-                                     * object since it makes no sense in 1d.
-                                     */
-    template <class VECTOR>
-    static void project (const Mapping<1>         &mapping,
-                        const DoFHandler<1>      &dof,
-                        const ConstraintMatrix   &constraints,
-                        const Quadrature<1>      &quadrature,
-                        const Function<1>        &function,
-                        VECTOR                   &vec,
-                        const bool                enforce_zero_boundary = false,
-                        const Quadrature<0>      &q_boundary = *invalid_face_quadrature,
+                        const Quadrature<dim-1>  &q_boundary = (dim > 1 ?
+                                                                 QGauss<dim-1>(2) :
+                                                                 Quadrature<dim-1>(0)),
                         const bool                project_to_boundary_first = false);
 
                                     /**
@@ -541,29 +529,11 @@ class VectorTools
                         const Function<dim>      &function,
                         VECTOR                   &vec,
                         const bool                enforce_zero_boundary = false,
-                        const Quadrature<dim-1>  &q_boundary = QGauss<dim-1>(2),
+                        const Quadrature<dim-1>  &q_boundary = (dim > 1 ?
+                                                                 QGauss<dim-1>(2) :
+                                                                 Quadrature<dim-1>(0)),
                         const bool                project_to_boundary_first = false);
 
-                                    /**
-                                     * Declaration of specialization
-                                     * of the previous function for
-                                     * 1d. At present, it is not
-                                     * implemented.
-                                     *
-                                     * The default value of the boundary
-                                     * quadrature formula is an invalid
-                                     * object since it makes no sense in 1d.
-                                     */
-    template <class VECTOR>
-    static void project (const DoFHandler<1>      &dof,
-                        const ConstraintMatrix   &constraints,
-                        const Quadrature<1>      &quadrature,
-                        const Function<1>        &function,
-                        VECTOR                   &vec,
-                        const bool                enforce_zero_boundary = false,
-                        const Quadrature<0>      &q_boundary = *invalid_face_quadrature,
-                        const bool                project_to_boundary_first = false);
-    
                                     /**
                                      * Create a right hand side
                                      * vector. Prior content of the
@@ -1243,14 +1213,6 @@ class VectorTools
                                       * Exception
                                       */
     DeclException0 (ExcNoComponentSelected);
-
-  private:
-                                    /**
-                                     * Null pointer used to
-                                     * denote invalid face
-                                     * quadrature formulas in 1d.
-                                     */
-    static const Quadrature<0> * const invalid_face_quadrature;
 };
 
 
index 98760a90c980c43549101e65ee4ce3d1f88f2627..e956d7248ea7d5ffb931c3331ff1d0cec4b8b32e 100644 (file)
@@ -276,79 +276,41 @@ VectorTools::interpolate (const DoFHandler<dim>           &dof_1,
 }
 
 
-#if deal_II_dimension == 1
-
-template <class VECTOR>
-void VectorTools::project (const Mapping<1>       &,
-                          const DoFHandler<1>    &,
-                          const ConstraintMatrix &,
-                          const Quadrature<1>    &,
-                          const Function<1>      &,
-                          VECTOR                 &,
-                          const bool              ,
-                          const Quadrature<0>    &,
-                          const bool              )
-{
-                                  // this function should easily be
-                                  // implemented using the template
-                                  // below. However some changes have
-                                  // to be made since faces don't
-                                  // exist in 1D. Maybe integrate the
-                                  // creation of zero boundary values
-                                  // into the project_boundary_values
-                                  // function?
-  Assert (false, ExcNotImplemented());
-}
-
-
-template <class VECTOR>
-void VectorTools::project (const DoFHandler<1>    &dof_handler,
-                          const ConstraintMatrix &constraints,
-                          const Quadrature<1>    &quadrature,
-                          const Function<1>      &function,
-                          VECTOR                 &vec_result,
-                          const bool              enforce_zero_boundary,
-                          const Quadrature<0>    &q_boundary,
-                          const bool              project_to_boundary_first)
-{
-  Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping"));
-  static const MappingQ1<1> mapping;
-  project (mapping, dof_handler, constraints, quadrature, function, vec_result,
-          enforce_zero_boundary, q_boundary, project_to_boundary_first);
-}
-
-
-#endif
-
-
-template <int dim, class VECTOR>
-void VectorTools::project (const Mapping<dim>       &mapping,
-                          const DoFHandler<dim>    &dof,
-                          const ConstraintMatrix   &constraints,
-                          const Quadrature<dim>    &quadrature,
-                          const Function<dim>      &function,
-                          VECTOR                   &vec_result,
-                          const bool                enforce_zero_boundary,
-                          const Quadrature<dim-1>  &q_boundary,
-                          const bool                project_to_boundary_first)
+namespace internal
 {
-  Assert (dof.get_fe().n_components() == function.n_components,
-         ExcInvalidFE());
+  namespace VectorTools
+  {
+#if deal_II_dimension == 1
 
-  Assert (vec_result.size() == dof.n_dofs(),
-          ExcDimensionMismatch (vec_result.size(), dof.n_dofs()));
-  
-  const FiniteElement<dim> &fe = dof.get_fe();
+    void
+    interpolate_zero_boundary_values (const ::DoFHandler<1>         &dof_handler,
+                                      std::map<unsigned int,double> &boundary_values)
+    {
+                                       // we only need to find the
+                                       // left-most and right-most
+                                       // vertex and query its vertex
+                                       // dof indices. that's easy :-)
+      for (unsigned direction=0; direction<2; ++direction)
+        {
+          ::DoFHandler<1>::cell_iterator
+              cell = dof_handler.begin(0);
+          while (!cell->at_boundary(direction))
+            cell = cell->neighbor(direction);
 
-                                  // make up boundary values
-  std::map<unsigned int,double> boundary_values;
+          for (unsigned int i=0; i<dof_handler.get_fe().dofs_per_vertex; ++i)
+            boundary_values[cell->vertex_dof_index (direction, i)] = 0.;
+        }
+    }
 
-  if (enforce_zero_boundary == true) 
-                                    // no need to project boundary
-                                    // values, but enforce
-                                    // homogeneous boundary values
-                                    // anyway
+#else
+    
+    template <int dim>
+    void
+    interpolate_zero_boundary_values (const ::DoFHandler<dim>       &dof_handler,
+                                      std::map<unsigned int,double> &boundary_values)
     {
+      const FiniteElement<dim> &fe = dof_handler.get_fe();
+
                                       // loop over all boundary faces
                                       // to get all dof indices of
                                       // dofs on the boundary. note
@@ -370,8 +332,9 @@ void VectorTools::project (const Mapping<dim>       &mapping,
                                       // that is actually wholly on
                                       // the boundary, not only by
                                       // one line or one vertex
-      typename DoFHandler<dim>::active_face_iterator face = dof.begin_active_face(),
-                                                    endf = dof.end_face();
+      typename ::DoFHandler<dim>::active_face_iterator
+        face = dof_handler.begin_active_face(),
+        endf = dof_handler.end_face();
       std::vector<unsigned int> face_dof_indices (fe.dofs_per_face);
       for (; face!=endf; ++face)
        if (face->at_boundary())
@@ -385,24 +348,60 @@ void VectorTools::project (const Mapping<dim>       &mapping,
                                               // vector valued elements here,
                                               // since we set all components
              boundary_values[face_dof_indices[i]] = 0.;
-         };
+         }
     }
+    
+#endif
+  }
+}
+
+
+
+template <int dim, class VECTOR>
+void VectorTools::project (const Mapping<dim>       &mapping,
+                          const DoFHandler<dim>    &dof,
+                          const ConstraintMatrix   &constraints,
+                          const Quadrature<dim>    &quadrature,
+                          const Function<dim>      &function,
+                          VECTOR                   &vec_result,
+                          const bool                enforce_zero_boundary,
+                          const Quadrature<dim-1>  &q_boundary,
+                          const bool                project_to_boundary_first)
+{
+  Assert (dof.get_fe().n_components() == function.n_components,
+         ExcInvalidFE());
+
+  Assert (vec_result.size() == dof.n_dofs(),
+          ExcDimensionMismatch (vec_result.size(), dof.n_dofs()));
+  
+                                  // make up boundary values
+  std::map<unsigned int,double> boundary_values;
+  
+  if (enforce_zero_boundary == true) 
+                                    // no need to project boundary
+                                    // values, but enforce
+                                    // homogeneous boundary values
+                                    // anyway
+    internal::VectorTools::
+      interpolate_zero_boundary_values (dof, boundary_values);
+  
   else
                                     // no homogeneous boundary values
     if (project_to_boundary_first == true)
                                       // boundary projection required
       {
-                                        // set up a list of boundary functions for
-                                        // the different boundary parts. We want the
-                                        // @p{function} to hold on all parts of the
-                                        // boundary
+                                        // set up a list of boundary
+                                        // functions for the
+                                        // different boundary
+                                        // parts. We want the
+                                        // @p{function} to hold on
+                                        // all parts of the boundary
        typename FunctionMap<dim>::type boundary_functions;
        for (unsigned char c=0; c<255; ++c)
          boundary_functions[c] = &function;
        project_boundary_values (dof, boundary_functions, q_boundary,
                                 boundary_values);
-      };
-
+      }
 
                                   // set up mass matrix and right hand side
   Vector<double> vec (dof.n_dofs());
index a4426aee0d03cf3c213ef333d9765068d104d36f..4b25b8524b82fab65d97524f62e5ec0ebf05b12a 100644 (file)
@@ -786,6 +786,13 @@ inconvenience this causes.
 <h3>deal.II</h3>
 
 <ol>
+  <li> <p>
+       Extended: The <code class="member">VectorTools::project</code> functions
+       are now also implemented for 1d.
+       <br> 
+       (WB 2006/08/08)
+       </p>
+
   <li> <p>
        Extended: <code 
                 class="class">DerivativeApproximation</code> now offers access to the

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.