]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Make the VectorTools::interpolate_boundary_values function available for the codim...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 2 Nov 2010 22:31:58 +0000 (22:31 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 2 Nov 2010 22:31:58 +0000 (22:31 +0000)
git-svn-id: https://svn.dealii.org/trunk@22588 0785d39b-7218-0410-832d-ea1e28bc413d

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

index ea5e92646e51d32d578890d52877c20a0dd082f0..c6e2146b5811ad594f4add4da9252971c7e23cf0 100644 (file)
@@ -332,6 +332,15 @@ through DoFHandler::get_tria() and DoFHandler::get_fe(), respectively.
 <h3>deal.II</h3>
 
 <ol>
+  <li><p>New: The class hp::FEFaceValues and hp::FESubfaceValues were not
+  previously available for meshes that were embedded in a higher
+  dimensional space (i.e. if the codimension was greater than 1). This is
+  now fixed. As a consequence, the VectorTools::interpolate_boundary_values
+  function is now also available for such meshes.
+  <br>
+  (WB, 2010/11/02)
+  </p></li>
+
   <li><p>Fixed: The FEValuesExtractors::Vector class did not work when the dimension
   of the domain was not equal to the dimension of the space in which it is
   embedded. This is now fixed.
index e0494237be9a9b360329dc8ebcee6061d51344c6..b27f75adf4f685d47da4523caeb7c70c630eaf01 100644 (file)
@@ -1353,6 +1353,7 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
                              const std::vector<bool>       &component_mask_)
 {
   const unsigned int dim=DH::dimension;
+  const unsigned int spacedim=DH::space_dimension;
 
   Assert ((component_mask_.size() == 0) ||
          (component_mask_.size() == dof.get_fe().n_components()),
@@ -1373,7 +1374,7 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
   const unsigned int        n_components = DoFTools::n_components(dof);
   const bool                fe_is_system = (n_components != 1);
 
-  for (typename FunctionMap<DH::space_dimension>::type::const_iterator i=function_map.begin();
+  for (typename FunctionMap<spacedim>::type::const_iterator i=function_map.begin();
        i!=function_map.end(); ++i)
     Assert (n_components == i->second->n_components,
            ExcDimensionMismatch(n_components, i->second->n_components));
@@ -1391,7 +1392,7 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
   std::vector<unsigned int> face_dofs;
   face_dofs.reserve (DoFTools::max_dofs_per_face(dof));
 
-  std::vector<Point<DH::space_dimension> >  dof_locations;
+  std::vector<Point<spacedim> >  dof_locations;
   dof_locations.reserve (DoFTools::max_dofs_per_face(dof));
 
                                   // array to store the values of
@@ -1411,11 +1412,11 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
                                   // the interpolation points of all
                                   // finite elements that may ever be
                                   // in use
-  hp::FECollection<dim> finite_elements (dof.get_fe());
+  hp::FECollection<dim,spacedim> finite_elements (dof.get_fe());
   hp::QCollection<dim-1>  q_collection;
   for (unsigned int f=0; f<finite_elements.size(); ++f)
     {
-      const FiniteElement<dim> &fe = finite_elements[f];
+      const FiniteElement<dim,spacedim> &fe = finite_elements[f];
 
                                       // generate a quadrature rule
                                       // on the face from the unit
@@ -1483,9 +1484,9 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
                                   // hp::FEFaceValues object that we
                                   // can use to evaluate the boundary
                                   // values at
-  hp::MappingCollection<dim> mapping_collection (mapping);
-  hp::FEFaceValues<dim> x_fe_values (mapping_collection, finite_elements, q_collection,
-                                    update_quadrature_points);
+  hp::MappingCollection<dim,spacedim> mapping_collection (mapping);
+  hp::FEFaceValues<dim,spacedim> x_fe_values (mapping_collection, finite_elements, q_collection,
+                                             update_quadrature_points);
 
   typename DH::active_cell_iterator cell = dof.begin_active(),
                                    endc = dof.end();
@@ -1494,7 +1495,7 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
       for (unsigned int face_no = 0; face_no < GeometryInfo<dim>::faces_per_cell;
           ++face_no)
        {
-         const FiniteElement<dim,DH::space_dimension> &fe = cell->get_fe();
+         const FiniteElement<dim,spacedim> &fe = cell->get_fe();
 
                                           // we can presently deal only with
                                           // primitive elements for boundary
@@ -1526,14 +1527,14 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
            {
                                               // face is of the right component
              x_fe_values.reinit(cell, face_no);
-             const FEFaceValues<dim> &fe_values = x_fe_values.get_present_fe_values();
+             const FEFaceValues<dim,spacedim> &fe_values = x_fe_values.get_present_fe_values();
 
                                               // get indices, physical location and
                                               // boundary values of dofs on this
                                               // face
              face_dofs.resize (fe.dofs_per_face);
              face->get_dof_indices (face_dofs, cell->active_fe_index());
-             const std::vector<Point<DH::space_dimension> > &dof_locations
+             const std::vector<Point<spacedim> > &dof_locations
                = fe_values.get_quadrature_points ();
 
              if (fe_is_system)
@@ -1997,7 +1998,7 @@ VectorTools::project_boundary_values (const Mapping<dim, spacedim>       &mappin
     && ! excluded_dofs[dof_to_boundary_mapping[i]])
       {
        Assert(numbers::is_finite(boundary_projection(dof_to_boundary_mapping[i])), ExcNumberNotFinite());
-       
+
                                       // this dof is on one of the
                                       // interesting boundary parts
                                       //
index c0bf0b59f7de5c7fbe59a24dfe02da4cd4c91f53..c30a9f4c51c223d5b76d8590c133443fa4496031 100644 (file)
@@ -206,6 +206,101 @@ void VectorTools::interpolate_boundary_values (
   ConstraintMatrix                    &,
   const std::vector<bool>    &);
 
+#if deal_II_dimension < 3
+template
+void VectorTools::interpolate_boundary_values (
+  const DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const unsigned char,
+  const Function<deal_II_dimension+1>   &,
+  std::map<unsigned int,double>       &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const hp::DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const unsigned char,
+  const Function<deal_II_dimension+1>   &,
+  std::map<unsigned int,double>       &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const MGDoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const unsigned char,
+  const Function<deal_II_dimension+1>   &,
+  std::map<unsigned int,double>       &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const unsigned char,
+  const Function<deal_II_dimension+1>   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const hp::DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const unsigned char,
+  const Function<deal_II_dimension+1>   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const MGDoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const unsigned char,
+  const Function<deal_II_dimension+1>   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const Mapping<deal_II_dimension,deal_II_dimension+1>    &,
+  const DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const FunctionMap<deal_II_dimension+1>::type   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const Mapping<deal_II_dimension,deal_II_dimension+1>    &,
+  const hp::DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const FunctionMap<deal_II_dimension+1>::type   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const Mapping<deal_II_dimension,deal_II_dimension+1>    &,
+  const MGDoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const FunctionMap<deal_II_dimension+1>::type   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const FunctionMap<deal_II_dimension+1>::type   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const hp::DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const FunctionMap<deal_II_dimension+1>::type   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const MGDoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const FunctionMap<deal_II_dimension+1>::type   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+#endif
+
 template
 void VectorTools::project_boundary_values<deal_II_dimension>
 (const Mapping<deal_II_dimension>     &,
@@ -266,22 +361,6 @@ VectorTools::compute_no_normal_flux_constraints (const DoFHandler<deal_II_dimens
 
 
 
-// // Due to introducing the DoFHandler as a template parameter,
-// // the following instantiations are required in 1d
-// #if deal_II_dimension == 1
-// template
-// void VectorTools::interpolate_boundary_values<deal_II_dimension>
-// (const Mapping<1>         &,
-//  const DoFHandler<1>      &,
-//  const unsigned char,
-//  const Function<1>        &,
-//  std::map<unsigned int,double> &,
-//  const std::vector<bool>       &);
-// #endif
-
-// the following two functions are not derived from a template in 1d
-// and thus need no explicit instantiation
-//#if deal_II_dimension > 1
 template
 void VectorTools::interpolate_boundary_values
 (const Mapping<deal_II_dimension>    &,
@@ -299,7 +378,24 @@ void VectorTools::interpolate_boundary_values
  std::map<unsigned int,double>       &,
  const std::vector<bool>    &);
 
-//#endif
+#if deal_II_dimension < 3
+template
+void VectorTools::interpolate_boundary_values
+(const Mapping<deal_II_dimension,deal_II_dimension+1>    &,
+ const DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+ const FunctionMap<deal_II_dimension+1>::type &,
+ std::map<unsigned int,double>       &,
+ const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values
+(const Mapping<deal_II_dimension,deal_II_dimension+1>    &,
+ const DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+ const unsigned char,
+ const Function<deal_II_dimension+1>   &,
+ std::map<unsigned int,double>       &,
+ const std::vector<bool>    &);
+#endif
 
 
 template

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.