]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement VectorTools::interpolate_boundary_values for FunctionMap argument.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 23 Apr 2001 15:24:10 +0000 (15:24 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 23 Apr 2001 15:24:10 +0000 (15:24 +0000)
git-svn-id: https://svn.dealii.org/trunk@4464 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/vectors.h
deal.II/deal.II/source/numerics/vectors.cc
deal.II/doc/news/2001/c-3-1.html

index 24068ac11166046e26c017ccd6e58ffe86109041..0c556f34ad64fd5ed8e533b47ce7ea98f481d264 100644 (file)
@@ -439,24 +439,29 @@ class VectorTools
                                        const Function<dim>   &rhs,
                                        Vector<double>        &rhs_vector);
 
-//TODO:[WB] Update interpolate_boundary_values for use of FunctionMap.
-//  keep both, the functions using and the functions not using a FunctionMap.
                                     /**
-                                     * Prepare Dirichlet boundary conditions.
-                                     * Make up the list of nodes subject
-                                     * to Dirichlet boundary conditions
-                                     * and the values to be
-                                     * assigned to them, by interpolation around
-                                     * the boundary. If the
-                                     * @p{boundary_values} contained values
-                                     * before, the new ones are added, or
-                                     * the old one overwritten if a node
-                                     * of the boundary part to be projected
-                                     * on already was in the variable.
+                                     * Prepare Dirichlet boundary
+                                     * conditions.  Make up the list
+                                     * of degrees of freedom subject
+                                     * to Dirichlet boundary
+                                     * conditions and the values to
+                                     * be assigned to them, by
+                                     * interpolation around the
+                                     * boundary. If the
+                                     * @p{boundary_values} contained
+                                     * values before, the new ones
+                                     * are added, or the old one
+                                     * overwritten if a node of the
+                                     * boundary part to be projected
+                                     * on already was in the
+                                     * variable.
                                      *
-                                     * The parameter @p{boundary_component} corresponds
-                                     * to the number @p{boundary_indicator} of the face.
-                                     * 255 is an illegal value, since it is reserved
+                                     * The parameter
+                                     * @p{boundary_component}
+                                     * corresponds to the number
+                                     * @p{boundary_indicator} of the
+                                     * face.  255 is an illegal
+                                     * value, since it is reserved
                                      * for interior faces.
                                      *
                                      * The flags in the last
@@ -467,13 +472,13 @@ class VectorTools
                                      * specified by the default value
                                      * (i.e. an empty array), all
                                      * components are
-                                     * interpolated. If is different
-                                     * from the default value, it is
-                                     * assumed that the number of
-                                     * entries equals the number of
-                                     * components in the boundary
-                                     * functions and the finite
-                                     * element.
+                                     * interpolated. If it is
+                                     * different from the default
+                                     * value, it is assumed that the
+                                     * number of entries equals the
+                                     * number of components in the
+                                     * boundary functions and the
+                                     * finite element.
                                      *
                                      * It is assumed that the number
                                      * of components of the function
@@ -485,15 +490,35 @@ class VectorTools
                                      * information.
                                      */
     template <int dim>
+    static void interpolate_boundary_values (const Mapping<dim>            &mapping,
+                                            const DoFHandler<dim>         &dof,
+                                            const typename FunctionMap<dim>::type &function_map,
+                                            std::map<unsigned int,double> &boundary_values,
+                                            const std::vector<bool>       &component_mask = std::vector<bool>());
+
+                                    /**
+                                     * Same function as above, but
+                                     * taking only one pair of
+                                     * boundary indicator and
+                                     * corresponding boundary
+                                     * function. Calls the other
+                                     * function with remapped
+                                     * arguments.
+                                     *
+                                     * This function is there mainly
+                                     * for backward compatibility.
+                                     */
+    template <int dim>
     static void interpolate_boundary_values (const Mapping<dim>            &mapping,
                                             const DoFHandler<dim>         &dof,
                                             const unsigned char            boundary_component,
                                             const Function<dim>           &boundary_function,
                                             std::map<unsigned int,double> &boundary_values,
                                             const std::vector<bool>       &component_mask = std::vector<bool>());
-
+    
                                     /**
-                                     * Calls the @p{interpolate_boundary_values}
+                                     * Calls the other
+                                     * @p{interpolate_boundary_values}
                                      * function, see above, with
                                      * @p{mapping=MappingQ1<dim>()}.
                                      */
@@ -504,6 +529,19 @@ class VectorTools
                                             std::map<unsigned int,double> &boundary_values,
                                             const std::vector<bool>       &component_mask = std::vector<bool>());
 
+                                    /**
+                                     * Calls the other
+                                     * @p{interpolate_boundary_values}
+                                     * function, see above, with
+                                     * @p{mapping=MappingQ1<dim>()}.
+                                     */
+    template <int dim>
+    static void interpolate_boundary_values (const DoFHandler<dim>         &dof,
+                                            const typename FunctionMap<dim>::type &function_map,
+                                            std::map<unsigned int,double> &boundary_values,
+                                            const std::vector<bool>       &component_mask = std::vector<bool>());
+
+    
 //TODO:[WB] Update project_boundary_values for more components
                                     /**
                                      * Project @p{function} to the boundary
index 27aa6fa8e8c21a6f74a429be049d8ddf614281fe..b21614cd81285382f92a683c88dd8524f6c886b9 100644 (file)
@@ -572,7 +572,7 @@ VectorTools::interpolate_boundary_values (const Mapping<1>         &,
   const std::vector<bool> component_mask ((component_mask_.size() == 0) ?
                                          std::vector<bool> (fe.n_components(), true) :
                                          component_mask_);
-  Assert (count(component_mask.begin(), component_mask.end(), true) > 0,
+  Assert (std::count(component_mask.begin(), component_mask.end(), true) > 0,
          ExcComponentMismatch());
   
                                   // check whether boundary values at
@@ -623,6 +623,22 @@ VectorTools::interpolate_boundary_values (const Mapping<1>         &,
 };
 
 
+
+template <>
+void
+VectorTools::interpolate_boundary_values (const Mapping<1>              &mapping,
+                                         const DoFHandler<1>           &dof,
+                                         const FunctionMap<1>::type    &function_map,
+                                         std::map<unsigned int,double> &boundary_values,
+                                         const std::vector<bool>       &component_mask)
+{
+  for (FunctionMap<1>::type::const_iterator i=function_map.begin();
+       i!=function_map.end(); ++i)
+    interpolate_boundary_values (mapping, dof, i->first, *i->second,
+                                boundary_values, component_mask);
+};
+
+
 #endif
 
 
@@ -630,19 +646,24 @@ template <int dim>
 void
 VectorTools::interpolate_boundary_values (const Mapping<dim>            &mapping,
                                          const DoFHandler<dim>         &dof,
-                                         const unsigned char            boundary_component,
-                                         const Function<dim>           &boundary_function,
+                                         const typename FunctionMap<dim>::type &function_map,
                                          std::map<unsigned int,double> &boundary_values,
                                          const std::vector<bool>       &component_mask_)
 {
-  Assert (boundary_component != 255,
+                                  // if for whatever we were passed
+                                  // an empty map, return immediately
+  if (function_map.size() == 0)
+    return;
+  
+  Assert (function_map.find(255) == function_map.end(),
          ExcInvalidBoundaryIndicator());
 
   const FiniteElement<dim> &fe           = dof.get_fe();
   const unsigned int        n_components = fe.n_components();
   const bool                fe_is_system = (n_components != 1);
   
-  Assert (n_components == boundary_function.n_components,
+  Assert ((function_map.size() == 0) ||
+         (n_components == function_map.begin()->second->n_components),
          ExcInvalidFE());
 
                                   // set the component mask to either
@@ -651,7 +672,7 @@ VectorTools::interpolate_boundary_values (const Mapping<dim>            &mapping
   const std::vector<bool> component_mask ((component_mask_.size() == 0) ?
                                          std::vector<bool> (fe.n_components(), true) :
                                          component_mask_);
-  Assert (count(component_mask.begin(), component_mask.end(), true) > 0,
+  Assert (std::count(component_mask.begin(), component_mask.end(), true) > 0,
          ExcComponentMismatch());
 
                                   // field to store the indices
@@ -686,14 +707,15 @@ VectorTools::interpolate_boundary_values (const Mapping<dim>            &mapping
   Quadrature<dim-1> aux_quad (unit_support_points, dummy_weights);
   FEFaceValues<dim> fe_values (mapping, fe, aux_quad, update_q_points);
 
-  DoFHandler<dim>::active_cell_iterator cell = dof.begin_active(),
-                                       endc = dof.end();
+  typename DoFHandler<dim>::active_cell_iterator cell = dof.begin_active(),
+                                                endc = dof.end();
   for (; cell!=endc; ++cell)
     for (unsigned int face_no = 0; face_no < GeometryInfo<dim>::faces_per_cell;
         ++face_no)
       {
-       DoFHandler<dim>::face_iterator face = cell->face(face_no);
-       if (boundary_component == face->boundary_indicator()) 
+       typename DoFHandler<dim>::face_iterator face = cell->face(face_no);
+       const unsigned char boundary_component = face->boundary_indicator();
+       if (function_map.find(boundary_component) != function_map.end()) 
                                           // face is of the right component
          {
                                             // get indices, physical location and
@@ -705,10 +727,10 @@ VectorTools::interpolate_boundary_values (const Mapping<dim>            &mapping
            
            if (fe_is_system)
              {
-               boundary_function.vector_value_list (dof_locations, dof_values_system);
+               function_map.find(boundary_component)->second->vector_value_list (dof_locations,
+                                                                                 dof_values_system);
                
                                                 // enter into list
-               
                for (unsigned int i=0; i<face_dofs.size(); ++i)
                  if (component_mask[fe.face_system_to_component_index(i).first])
                    boundary_values[face_dofs[i]]
@@ -720,9 +742,9 @@ VectorTools::interpolate_boundary_values (const Mapping<dim>            &mapping
              {
                                                 // get only the one component that
                                                 // this function has
-               boundary_function.value_list (dof_locations,
-                                             dof_values_scalar,
-                                             0);
+               function_map.find(boundary_component)->second->value_list (dof_locations,
+                                                                          dof_values_scalar,
+                                                                          0);
                
                                                 // enter into list
                
@@ -733,6 +755,24 @@ VectorTools::interpolate_boundary_values (const Mapping<dim>            &mapping
       }
 }
 
+
+
+template <int dim>
+void
+VectorTools::interpolate_boundary_values (const Mapping<dim>            &mapping,
+                                         const DoFHandler<dim>         &dof,
+                                         const unsigned char            boundary_component,
+                                         const Function<dim>           &boundary_function,
+                                         std::map<unsigned int,double> &boundary_values,
+                                         const std::vector<bool>       &component_mask)
+{
+  typename FunctionMap<dim>::type function_map;
+  function_map[boundary_component] = &boundary_function;
+  interpolate_boundary_values (mapping, dof, function_map, boundary_values,
+                              component_mask);
+};
+
+
   
 template <int dim>
 void
@@ -750,6 +790,21 @@ VectorTools::interpolate_boundary_values (const DoFHandler<dim>         &dof,
 
 
 
+template <int dim>
+void
+VectorTools::interpolate_boundary_values (const DoFHandler<dim>         &dof,
+                                         const typename FunctionMap<dim>::type &function_map,
+                                         std::map<unsigned int,double> &boundary_values,
+                                         const std::vector<bool>       &component_mask)
+{
+  Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping"));
+  static const MappingQ1<dim> mapping;
+  interpolate_boundary_values(mapping, dof, function_map,
+                             boundary_values, component_mask);
+}
+
+
+
 template <int dim>
 void
 VectorTools::project_boundary_values (const Mapping<dim>       &mapping,
index 6fefb46dc0f59efecad5079eff05b35721d5544c..169aeb6aa5233b237219bb8068bce24b865dd198 100644 (file)
@@ -313,6 +313,16 @@ documentation, etc</a>.
 <h3>deal.II</h3>
 
 <ol>
+  <li> <p>
+       New: Implement the <code
+       class="member">VectorTools::interpolate_boundary_values</code>
+       a second time, this time taking a <code
+       class="class">FunctionMap</code> object as do all the other
+       functions of this type.
+       <br>
+       (WB 2001/04/23)
+       </p>
+
   <li> <p>
        Fixed: The vertex numbers written by the <code
        class="class">GridOut</code>::<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.